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On the Addressing Problem of 
Loop Switching 


By L. H. BRANDENBURG, B. GOPINATH, and R. P. KURSHAN 
(Manuscript received April 5, 1972) 


R. L. Graham and H. O. Pollak recently proposed an addressing scheme 
to J. R. Pierce’s Loop Switching Network-a network system which would 
exploit the one-way nature of much data transmission. Our paper provides 
a generalization of the scheme and proposes a new approach to certain 
unanswered questions. 


I. INTRODUCTION 


Recently, J. R. Pierce’ proposed a data communication system con- 
cept whereby subscribers are placed on interconnected “loops”. A 
message originating in one loop and destined for a recipient in another 
loop must find its way through the network of interconnections. Such 
a system is potentially efficient for data communication because it 
exploits the one-way nature of many such transmissions by eliminating 
the need to set up an entire line between two parties in the conventional 
manner. ; | 

The realization of this potential efficiency depends in part on the 
existence of an efficient routing scheme. Such a scheme is discussed by 
Graham and Pollak in Ref. 2. Briefly, they assigned to each loop a 
ternary address; a message destined for a certain loop is tagged with 
the address of its destination, and upon entering a loop, a computation 
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involving that loop address and the message destination address yields 
the distance of the message from its destination. The message then 
moves in such a way as to decrease this distance. In particular, it is 
assumed that there are n loops and each address is a sequence of N 
(fixed) 0’s, 1’s and d’s; let H be the n X N matrix the ith row H; of 
which is the address of the ith loop. The addresses are chosen in such 
a way that if D,; is the distance between the 7th and jth loop (by D;; 
is meant the fewest number of loops which it is necessary to cross to 
to get from the ith to the jth), then 


where a@ b = ‘b if the set {a, b} = {0, nh 
0 otherwise 


N 
D;; = Ss (Hi, ® Air), 
k=1 


(this is an “extended”? Hamming distance between H, and H;). 

Alternatively, the problem may be considered in graph-theoretic 
terms: associate with the loop network a graph the vertices of which 
correspond to the loops and the edges of which correspond to direct 
connections between two loops. Then D;,; is the minimal number of 
edges in a path between the 7th and jth vertices. 

Central to the efficiency of the above scheme was the possibility of 
keeping N small; if G is a connected graph on n vertices, then Graham 
and Pollak showed that an address matrix H with corresponding N—call 
it N,—could always be found to satisfy Ny, S (n — 1) X diam G (diam G 
= max,,; D;;). Unfortunately, (n — 1) X diam G can be unpleasantly 
large. However, they conjectured that, in fact, for any connected 
graph G, one could always find an H matrix such that Ny S n — 1. 
This was proved in the case of complete graphs, trees and cycles. 
Furthermore, defining N(G@) = min, N,,, they showed that, N(G) = 
n — 1 for G a cycle on an odd number of vertices, a complete graph or 
a tree; and N(G) = n/2 for G a cycle on an even number of vertices. 
In any case, it was shown that N(@) 2 max {n*, n-} where n*(n’) is 
the number of positive (negative) eigenvalues of the distance matrix D 
with elements D;; defined above. 

For a graph G' (which we will always assume to be connected) referred 
to implicitly or explicitly, let us reserve the letter D to denote the ° 
distance matrix of G, and n to denote the number of vertices of G. 
Then the conjecture of Ref. 2 becomes 


N(@ <n-1. (1) 


Incidentally, this conjecture is equivalent to the following conjecture: 
that any graph G on n nodes can be embedded in a “factor” F of an 
(n — 1)-cube in such a way that the distance matrix for G is a submatrix 
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of that for F; by a “factor” of a cube C is meant a graph F defined as 
follows: Let F = {f,, ---, fz} be a set of pairwise disjoint faces (of 
various dimensions) of C; let V be the set of vertices of C appearing 
in no f;. Then the set of vertices of Fis V U & (i.e., the vertices V of C 
along with the faces of 5). There is an edge between two vertices 2, 
ye VU Sil: 


(i) x, ye V and z is connected to y by an edge in C; 
(il) x, y e ¥ and some vertex of (the face) x is connected by an edge 
to some vertex of (the face) y; 
(iii) 2 e V, ye & and some vertex of (the face) y is connected by an 
edge to x. 


More descriptively, F is formed from C by shrinking each f; to a vertex. 

In this paper, we discuss the conjecture (1) by formulating the ad- 
dressing problem in matrix-algebra terms and then present an algebraic 
generalization of this that immediately yields an addressing scheme 
more efficient than that of Ref. 2. 

Specifically, we show that a graph G has an n X N address matrix 
H if and only if there exist n X N binary (0, 1)-valued matrices A 
and B such that 


D= ABt+ BAt (2) 


(where t means ‘‘transpose”’). Using this, we expand the class of graphs 
known to satisfy (1). Then we rewrite (2) as 

B* 

D = (A: Bil-- (3) 

At 
(the product of the two pairs of concatenated matrices). Note that 
the problem of finding n X N matrices A and B satisfying (2) [and (3)] 
can be generalized to the problem of finding n X 2N binary (0, 1)-valued 
matrices P and Q satisfying 


D = PQ. , (4) 


Clearly, any decomposition (2) yields the corresponding decomposition 
(4). On the other hand, any decomposition of the form (4) with P and Q 
(binary-valued) n X JZ matrices gives rise to a binary addressing of 
length L which, as in Ref. 2, will enable a message to compute its distance 
from its destination. In particular, the 7th vertex is tagged with Q; , 
the zth row of Q, while a message destined for the 7th vertex is tagged 
with P; , the 7th row of P. Then the message, at vertex 7, makes the 
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computation P,Q; = D;,; and determines its distance from its destina- 
tion. This computation involves summing JL products, each of the 
form a-b where a, be{0, 1}, as compared to summing N sums of the form 
a @ b defined above where a, b, e{0, 1, d}. Note that two sets of addresses 
are required: the set {Q;} is stored in the network, and members of the 
set {P;} are tagged onto messages. But users of the network only need 
a directory of the P,’s. 

Hence, assuming that the addresses of Ref. 2 are to be transmitted 
in binary, the decomposition (4) gives rise to no larger, and in certain 
cases, smaller addresses than those derived in Ref. 2. If L(G) is the 
minimum value of L over all such decompositions (4), then corresponding 
to the conjecture (1) [and a necessary consequence, if (1) is true] we 
conjecture that 


L(G) S 2(n — 1). (5) 
In any case, we show that 
L(@) s Dy max D,;; S n(diam G.) (6) 


This improves on the upper bound derived in Ref. 2, which after trans- 
lation of addresses into the form of (3), gives LD S 2(n — 1) diam G. 

Using the decomposition (4) which yields (6), we discuss a secondary 
coding of the addresses which reduces the number of bits actually 
transmitted; by using this decomposition rather than one with possibly 
smaller L, we achieve an addressing scheme the addresses of which are 
very easy to compute and which is hence amenable to continual up- 
dating of the loop network. 


II. EQUIVALENT FORMULATIONS OF THE GRAHAM—-POLLAK SCHEME 


Here we present three theorems: the first characterizes in various 
terms the property that a graph G admits an n X WN address matrix H, 
the second and third present some sufficient (but not necessary) con- 
ditions for (1) to hold. 

We will call a matrix binary-valued if its entires are either 0 or 1. 

Clearly, G admits an n X N A if and only if N = N(G). 


Theorem 1: Let D be the distance matrix of a graph G on n vertices. Then 
the following conditions are equivalent: 

(2) N 2 NG); 

(iz) There exist binary-valued n X N matrices A and B such that 
D = ABt + BA’; 
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(122) There exist n X N matrices C and E, C binary-valued and EF with 
entries 0, 1 such that C;; = | E;; | and 2D = CCt — EE*. 

(2v) There existn X n matrices F; ,1 = 1, -+- , N with each F; binary 
valued, symmetric, rank <2, and D = >.™_, F;. 


Proof: (2) = (iz): Since N 2 N(G) there exists for Gann X N address 
matrix H. By construction we know that D,,; is the extended Hamming 
distance between H,; and H, , which is in fact the number of positions 
where H; and H; differ without either being a d. Define for each H; a 
binary-valued N-vector (i.e., 1 X N matrix) A; by Ay, = lif Hy, = 1 
and A,;, = 0 otherwise, and define a binary-valued N-vector B; by 
B,, = 1if H,, = 0 and B;, = 0 otherwise. It is then clear that D,; is 
the sum of inner products (A; , B;) + (B;, A;), ie., Di; = A:(B;)* + 
B,(A,)*. We thus obtain D = ABt + BA'*, where A[B] is the n X N 
matrix whose 7th row is A,[B;]. 

(22) => (4): Suppose we have n X WN binary-valued matrices 
A and B such that D = ABt + BA’. For all 7, since D;; = 0, A;(B,)t = 
0 = B,(A,)* and A; and B; can hence have no 1’s in common positions. 
It should be clear from the converse case above that by defining H as 


i if A;; al | 


1 
H,; = 90 if B;; = ] ) 
li if Aj; = 5,; = ;l 
we obtain the desired address matrix. 

(17) => (vit): Given (iz), define matrices C and HasC = A+B 
and # = A — B. Since A and B have no common ones, C is binary- 
valued with zeros only in those positions where A and B have common 
zeros and one elsewhere. Moreover, # has values 0, +1 and has zeros 
only in those positions where A and B have common zeros and is -1 
elsewhere. The definition of C and E immediately leads to verification 
of the relation 

2D = CCt — EE*t. 


(22t) => (a2): Given (272), define A = 4(1 + EF) and B = 3(C — E). 
Clearly, A and B are binary-valued (with no ones in common positions) 
and (72) follows directly from (772). 

(i) = (wv): Let H’,7 = 1, --- N be the zth column of matrix 7 
and let F; be the n X n matrix of Hamming distances between the 
rows (actually sequences of length one) of H’. Then clearly 


N 
D = > F; 
i=1 
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and we must show that each Ff; has the required structure. Suppose 
that H* has n, 0’s, n2 1’s, and n; d’s. Reorder the elements of H’ so that 
all of the 0’s appear first, followed by all of the 1’s. Then, the rows and 
columns of F’; can be reordered (symmetrically) to yield 














ny No Ng 
[O---0 | 1---1 | 0-- 0 | 
| ae 
028504 eee || Sete 
Leeed OO) aiee8 
Ng = F, (7) 
ete ew hs Kee? 
es.) .eawd 0---0 
rl x 
Osa stees 0 eo 





Clearly, F; is binary-valued, symmetric, and rank two, (or rank zero 
if nm; or m2 = 0) and these properties are invariant with respect to 
symmetric permutation of rows and columns. 

(iv) => (i): Let D = >>™_, F; with each F, satisfying the condi- 
tions of the theorem. We show that by appropriate permutation of 
rows and columns, /’; can be partitioned as in (7). Note first that since 
diagonal elements of D are zero, the same must be true of each F,. 
Hence rank (7) S 2 implies either rank (F;) = 2 or F; = 0 since F; is 
symmetric. 


Now, let F be any binary-valued, symmetric, rank two, n X n 


matrix with zeros on the main diagonal. Let e,, --- e, represent the 
diagonal elements of F. (The value of each e; is 0). Two elements e; and 
e; will be defined as equivalent, e; ~ e; if F;; = 0. We now show that 


under the assumptions concerning F’, ~ is a true equivalence relation 
with respect to the e;’s that do not index a row or column of zeros in F'. 
Clearly, e; ~ e;,7 = 1, --- n, and e; ~ e; implies e; ~ e; . In order to 


t The equivalence class method of proof was suggested to the authors by H. 8. 
Witsenhausen. 
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prove transitivity, let us (temporarily) remove from F all rows and 
columns consisting entirely of zeros. Suppose now that e; ~ ¢; ,¢€; ~~ & , 
but that F;, = 1. Then F has the form 


























ce @ m 
[ 7] 
1 0;0/}1 y 
9 0;0)10 1 
k 1/0/;0 
m 0 ee: 0 
L 4 
Since rows and columns of zeros have been removed, F;,, = F',; = 1 


for some m. The unspecified values of F’,,, and F;, are denoted by x 
and y respectively. Clearly the 3 X 3 submatrix of F 


O01 y 
0 0 1 
1 x O 


consisting of columns j, k, m and rows 2, 7, m, is nonsingular regardless 
of how z and y are specified, thus contradicting the rank two condition. 
Hence, ~ is transitive, and F can be permuted in such a way that all 
diagonal blocks consist entirely of zeros, and all off-diagonal blocks 
consist entirely of ones. This partitioned matrix can then be bordered 
on the right and bottom with the previously removed rows and columns 
of zeros. The resulting matrix then resembles (7) except that we must 
discount the possibility of having more than two equivalence classes. 
However, it is obvious from the partitioned structure of F that the 
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number of equivalence classes equals rank (F) which is two by assump- 
tion. To each F’; we can thus associate a column vector H" consisting of 
0’s, 1’s and d’s. We place d’s in positions of H* corresponding to the 
elements e, that index rows and columns of 0’s in F;. We place 0’s 
in positions of H* corresponding to the elements e, that determine one 
of the equivalence classes, and we place 1’s in the remaining positions 
of H’. The matrix of Hamming distances associated with H* is thus 
equal to F;. If rank (F;) = 0 for some 7, then F; = O and the corre- 
sponding column H’ can be chosen to consist entirely of d’s. This com- 
pletes the proof. 

We now present a method for joining address matrices of graphs to 
form an address matrix of their union. Suppose a graph G contains 
(connected) subgraphs G, and G, such that G = G, U G, and G, 1 G, 
are exactly a vertex of G (i.e., G is found by connecting the two graphs 
G, and G, at some vertex). Then we will say that G is separable into 
G, and G.. 

Theorem 2: Suppose the graph G is separable into the subgraphs G, 
and G,. Then N(G) S N(G,) + N(G,). 

Proof: Let n; be the number of vertices of G; and let H; be ann; X 
N(G;,) address matrix for G; for = 1, 2. We may assume that G, 1 G, 
is the 2,th vertex of G, and the lst vertex of G,. Define H as ° 


N(G,) N(G2) 

— -, 

H,| Jz 
Ny 

H =\|—— 

A ee 


where J,(J2) is the (nx — 1) X N(Gi)((m — 1) X N(G,)) matrix 
each row of which is the n,th (1st) row of H,(H,). Clearly H is an ad- 
dress matrix for G with addresses of length N(G,) + N(G2). 





Corollary: If the graph G is separable into subgraphs each of which 
satisfies the conjecture (1), then G also satisfies (1): 1.e., N(G) S n — 1. 


Proof: With notation as above, if N(G,) S m — land N(G.) S nz —1 
then N(G) Sn+t(m—-1)-l=n-—1. 


Clearly, this corollary can be extended to any number of subgraphs. 
We will now present a class of graphs for which N(G) S n — 1. Call 
a linear tree (i.e., Di; = |i — 7 |) contained in a graph G a spine if each 
vertex of G is either a vertex of the linear tree or else distance one away 
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- from the linear tree. Clearly, any diameter two graph (i.e., the maximum 
distance between each pair of vertices is two) admits a spine; any 
maximal linear tree will do. It can also be shown that any diameter 
three graph admits a spine, while, on the other hand, there exists a 
14-point nonseparable diameter four graph admitting no spine. 

If L is a spine for a graph G, label the vertices of LZ consecutively: 
Vio, V20 y 30 ++ + and label each vertex v of G not on Lv,;,j = 1, 2, --: 
where 2 = min {k | d(v, v0) = 1}. 


Theorem 3: A sufficient condition for N(G) S n — 1 is that G admit 
a spine L satisfying 


EITHER (1) dQv;. , viz) 2 dio, Viz) — yx for ¢ = 1, 
k,e = 0, 1, --+ (when defined) where y, = 0 if k 
and1liafk € 0. 

OR (2) There are no six vertices x11, Y1, Yo, Lo, X3, %4, distinet 
except possibly y, = Yo, and increasing from left to right 
in the lexicographic order of the vertices’ double-indexing, 
with each x; a vertex of L and such that: 


2) d(te , a1) = Aye, %1) + 1 + y(yo), 
Ut) (te 5 Lo) = dee, Xo) + 1, 
11) U(r, , 21) = d(as, 21) + 1, 

wv) d(xs , Yr) dO 54) FP); 


where y(y:) = 0 7f y; is a vertex of L and 1 otherwise. 


2, POS od 
0 


| V 


[Note: (1) = (2), and there exist even weaker (but more complicated) 

sufficient conditions than (2); however, no graph 1s known to the authors 

which admits a spine and farls to satisfy (1) for some spine]. 

Corollary 1: For any diameter 3 graph G admitting a spine of 3 or fewer 
vertices, N(G) S n — I. 

Corollary 2: For any diameter 2 graph G admitting a spine L such that, 
either no four vertices of L form a 4-cycle in the original 
graph, or else L has at most 4 points, we have N(G) S n — 1; 
this includes all complete bipartite graphs. 


Proof: Both statements follow directly from condition (2) of the 
theorem. 


Proof of Theorem: Suppose we have a graph which admits a spine; 
double-indexing the vertices as above, the associated distance matrix D 
will have the form 
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Vio Yi Vie s** Va Vaeo *** U39 031032 *** Vag Vaan 
11 ie eee <2 | <3 | <4 
|| | | 
0 x2] | ! 
<3 | <4 | <5 
| | i] 
| | I 
ao) 3 ! 
rO4d hese g yd 2 | <3 
| | | | | | 
1 110 
i i i 
$3 i1 <2 | <3 <4 
Rl aati ! 
\ <2 
Dees Wee Maks cae Pak ome a Oe ee eet eee 
eae | | O|1 1 [eto <2 
[| —_} -——, —_|—______ | 
| | 1/0 S21 
<4 | <2 = a a <3 
a 
; <2 
1 1 0; 
ee Se Se We eet fe es eee ee ee Wine oe ee a er ee Se eee 
; i 1 tr OFLA 
[| | | 
3 | Pi fo a 
<5 | $3 <4 i 2 <3 ! 1 
; en 
1 I i] 





Conditions (1) and (2) both place certain restrictions on the extent of 
non-monotonicity of any column of D below the main diagonal; (1) says 
that the columns under each v,, shall be monotonically increasing and 
that the other columns shall contain decreases of at most one unit per 
square (roughly speaking). Condition (2) allows a greater degree of 
non-monotonicity but excludes an unmanageable combination of such. 
It is not difficult to see that (1) = (2). 

Now assume that (2) is satisfied. We will construct n X (n — 1) 
matrices A and B satisfying D = ABt + BA+*. Roughly speaking, one 
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can think of A as being superimposed on the n X (n — 1) matrix D’ 
remaining after removing from D its first column; we will use the 
notation [k, , k. , ---] to indicate that the designated spot in A shall 
be a 1 if the associated spot in D’ is either k, or k. or --- ; all non-1 
entries in A are 0. The matrix Bt will basically be ann — 1 KX n — 1 
identity matrix with a concatenated first column of zeros; 1’s will then 
be added above the existing diagonal of 1’s for the purpose of ‘“repeating”’ 
the associated column of A. 

First, consider the case in which (1) holds. In this case (it is easily 
checked) ABt can have the form 


1 1- 7 1 [2] ! [3] | [4] 
0 ci ! , 
11} 3} [4] | [5] 
2 S| es oe ee ere 
! raeet i [2] ! [3] 
eee 3 
2,3) | yop of fl hw 
Bre ! 
eee Je ee Sp eee 
! 0 1. T< 11 1 [2] 
! ‘olo 4 
! Bale od” hae 
[2, 3, 4] | (2,3) 1 | ee oe jae ©) 
! ja: 
a eee 2 os ee 
! 1 0 i 0 11 
0 0 ! 0 1 
[2, - hee [2, 3, 4] - [2, 3] | ay 
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| i] | 
OE Vere 1 IS 8520 51,952,353. 
I | | 
Alb Pi 
i 1, 2-1 
* | | | 
° cd) heey 3-32 | 
*: [oe { t 
- | 21 rs are 43 
0 0:1 
ee | ee Vics IN se ce oa a a Ee LS a ee ee es 
0 iO/1 1 111 
| 1 
~s a _ as 
1/0 a. ! 
| | { 
0) orgy 4 i beh 4 
‘ (jee | | 
Jie ee a | 
~ 0 iJ O11 ! 
0 \ J 0 11-1 11 
—— [=F Ie. (—_ 
| | | 
: ! op a ae 
| 1 { 
: 23,41} | 2,31 ie i 
; 2-1 = 
0 ! 1 1 0 11 
0 il 1 0 11 
| | 
0 | | 11 0 
1 | 
0 it 
| i} l 
[ if (3% 
l t t 
| | | 
where k,k, --- — 1 means that whenever k, or k, or --- appears in D, 


then J appears in ABt*. It is then clear that ABt + (AB*)t = D. 

The general case of (2) now follows in an analogous fashion; the 
matrix Bt will be as above, except that certain of the off-diagonal 1’s 
will be replaced by 0’s. The matrix A will be as above, except that the 
bracket notation will have to be written for each column above the 
diagonal, each row below the diagonal and the entries in the brackets 
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will have to be altered to meet the needs of the situation. It can be 
shown that given the conditions of (2), this can actually be done. 


Ill. EXTENSION TO THE PQt SCHEME 


As explained in the introduction, an alternative way of looking at 
the addressing problem is to determine the smallest L for which there 
are n X JL binary-valued matrices P and Q such that D = PQ*. This 
formulation also allows for a possibility where D is nonsymmetric. A 
nonsymmetric D can arise if the given graph is directed, or can arise 
if the elements of D are chosen so as to obtain a desirable routing 
scheme not strictly dependent on the actual distances between vertices. 
In these terms, the analog of the Graham and Pollak conjecture (1) is 
that L(G) < 2(n — 1) (where L(G) is the smallest ZL). We know this 
is true when V(G) S n — 1 since in any case (from (3) and (4)) 


L(G) S 2N(@). (8) 
In Theorem 4 we will show that L(G) S )>°; max; D;; < n diam G. 


On the other hand, we have the lower bounds 
L(G) 2 rank D (9) 
and 
L(G) 2 2 diam G. (10) 


The bound (9) follows since L(G) 2 rank P = rank D. To see (10), 
find 7 and j such that D;; = diam G. Then P;(Q;)* = diam G = P;(Q,;)*. 
Hence P; and Q; have diam G 1’s in common and the same is true for 
P; and Q; . But D;; = 0, so that in each spot that P; has a 1, Q; has a 
zero. Thus, Q; must contain at least diam G 1’s and diam G 0’s. 
Theorem 4: L(G) S >>; max; D;; S$ n diam G 


Proof: Let s; = max; D;; and r, = >.'_; s; with r, = 0. It is easy to 
see that the following construction for P and Q gives an address of 
length r, . For every k, the kth row of P, P, is zero everywhere except 
in positions r,-; + 1 to 7, , where it is one. For every 7, the jth row of Q, 
Q, is divided into cells, the kth cell defined as positions 7,_, + 1 tor, . 
In the kth cell is placed exactly D,; 1’s with the rest being equal zero. 
Obviously P,Q; will then be the number of 1’s in the kth cell of Q; ., 
i.e., D,; . Notice that this proof is not dependent on D being the distance 
matrix of a graph; the theorem is true for any nonnegative integer-valued 
matrix D. 

A comparison of (3) and (4) suggests that the inequality (8) can be 
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strict. This is illustrated in several of the following special cases where 
we show, In particular, that if G is a cycle with n odd or a complete 
graph, then L(G) = N(G) + 1, and for certain trees, L(G) S $N(Q). 


(2) Complete Graphs 


For a complete graph, all off-diagonal elements of D are 1’s. Let 
P = TI, the n-dimensional identity matrix, and let Q = D. Then ob- 
viously D = ID = PQ+*, and in light of (9), since rank D = n, this 
factorization 1s minimal. 

For a binary-valued matrix X, define X to be X with all the 1’s 
changed to 0’s and all the 0’s changed to 1’s. Let X* = X*t. Then the 
decomposition for the complete graph above has the special form 
D = PP*, This property has important implications which will be 
discussed later. 

(it) Trees 

In Ref. 2 (Theorem 3), it was shown that for any tree G, N(G) = 
n — 1, and in fact, the address matrix H could be chosen to be binary- 
valued-without d’s. Hence, by (3), for any tree G, L(G) S 2(n — 1). 
For G a linear tree, 1.e., D;; = |¢ — 7 |, in fact L(G) = 2(n — 1) by (10), 
and hence for an arbitrary graph G there is no general upper bound for 
L(G) lower than 2(n — 1). However, we shall show that for some trees 
L(G) < 2(n — 1). Also, since the address matrix for a tree can be 
chosen to have no d’s, it is clear from (3) and Theorem 1 that the dis- 
tance matrix for a tree can be written as D = PP*. 

Now let G be a star on five vertices. 








Z 3 
1 
4 
5 
Then 
fo 111 1] 
Po 2 BB 
D=/1 2 0 2 2 
eee ae a 
l1 2220] 
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Let 
100100 
OL 0-1-0 4 
P=/1 0100 1}, 
11001 0 
O00 a 6 
and 
0100 4] 
10101 
Fee eae a 
00110 
01100 
00011 


Then D = PQ* (note that we no longer have Qt = P*) and L(G) S 
6 < 8 = 2N(G). Conversely, we show L(G) = 6. To see this, suppose 
there exists P and Q of order 5 X 5 such that D = PQ*. Then for D, , 
the submatrix of D defined below we have 


02 2 2 
Dow le 2 2) pete eRe 
220 2 


22 2 0 


where R and S are the obvious submatrices of P and Q, respectively. 
Observe first that each row of R and S must contain at least two I’s. 
Also, at least one of the rows of R or S must contain at most two 1’s, 
since if this were not the case, the diagonal elements of D, would not 
vanish. Since the roles of R and S are interchangeable, we assume with- 
out loss of generality that the first two elements of the first row R, of R 
are 1’s and the rest are zeros. Now since 
RS; = 2, J = 2, 3, 4, 

we have 

Si = Se = 1, j = 2,3, 4; 
thus 


Ra = Rye = j = 2,3,4 
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because 
R,;S; = 0 for all j. 


Now, at least two of the numbers R,; ,7 = 3, 4, 5 must be 1, and hence, 
at most one of the numbers S.; , 7 = 3, 4, 5 can be 1, which implies 
R;S; = <1, which contradicts the condition that R,;S; = 2. Thus 
L(G) > 5. 

If G is a star with n = 4k + 1, then L(G) S 6k < 8k = 2N(G). 
This follows from the above example and a result analogous to Theorem 
2: viz., if G is separable into Gj; and G, then L(G) S L(G,) + L(G). 


(ziz) Cycles—n odd 

Suppose G is a cycle with distance matrix D and n odd. Let m = 
(n — 1)/2 and let L,, be the distance matrix for the linear tree on im 
vertices. Let A be the m X m triangular matrix with 1’s on and below 
the diagonal, and let U be the m X m matrix all the entries of which are 1: 





[1 Lb chee. Ze 
A= “= ee ee 
1 
L 1 1 1 
Then 
m 
m—1 
Dy, mU — L, + A* 
2 
1 





D={;m m—1::-2 1 0 1 2-+-m—1 mi: 


mU —L, +A Ly, 
m—1 


m 
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Define 





Then, by direct computation, we get 


AA*+A*A | Dy At| AA‘ + A‘A* 





PP* = » (A* + 0) 0 pp. 
AA*+AA+U)| O'A | AA‘ + AA* 


where >.’(>.°) means: sum the rows (columns) of, and I is the m X m 
identity matrix. It can be shown that PP* = D. Note that Pisn X n 
and hence L(G) S n; but D is nonsingular (see Ref. 2, Theorem 4) 
and so by (9), L(G) = n. 


(iv) Cycles—n even 

In Ref. 2 (Theorem 4), it was shown that for G, a cycle with n even, 
N(G) = n/2, and hence by (8), L(G) S n; but diam G = n/2 and by 
(10) L(G) = n, so L(G) = n. Furthermore, the minimal address matrix 
H for G can be chosen without d’s and hence we can write D = PP* 
with Pn X n. 


(v) Cube 


Let G be the vertex and edge structure of the m-dimensional cube; 
then n = 2”. It can be shown that N(G@) = m, and since diam G = m, 
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L(G) = 2m. Also, it is apparent that we can write D = PP* with 
Pn X 2m. 


(vt) Diameter 2 Graphs 


We assume that the graph has at least four vertices, since a graph on 

three or fewer vertices must be a complete graph or a tree, and we know 
how to address such graphs. There are now two cases to be considered. 
Case 1. Here, we assume that the graph is a star. We know N(G) = 
n — 1, so L(G) S 2(n — 1). 
Case 2. Here we assume that the graph is not a star. Then, since the 
graph is connected, there must be four vertices, say 0; , V2, Un-1 , and », , 
such that v, is adjacent to v2 and v,_,; 1s adjacent to v, . The distance 
matrix has the structure 


0 1 
1 0 


0 1 
1 0 
where the unspecified elements are 0’s, 1’s, and 2’s. Let D, represent 


the distance matrix of a complete graph on n vertices and let D = 
D, + D, . The remainder matrix D, has the form 


0 0 


0 0 


0 0 
l 0 0 


and the unspecified elements are 0’s and 1’s. The idea now is to take 
1’s from D, and put them into D, in such a way as to reduce the rank 
of D, . This can be done by taking ones from the first two and last two 
rows of D, and placing them into the first two and last two rows of D, . 
The result is D = Di + Dj where 
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1 1;1---1;0 0 
Lot dk ee Oe 


and D/ remains a 0, 1 valued matrix. Note that rank (D{) < n — 2. 
Let J,-2 and J, be the identity matrices of order n — 2 and n respec- 
tively. Then D = PQ*, where 


ar 


0:--0 1 
and 
—— SS 
Di ‘lin 





1 1/1---1/0 OJ 


Hence, L(G) S 2(n — 1). It is not yet known whether N(G) S n — 1 
for all diameter 2 graphs. Note that in the above construction P* #~ Qt 
in general. 


IV. SOME LINEAR ALGEBRA 


In the previous section we realized the representation (4) for several 
classes of graphs. It is possible to evaluate lower bounds for N(G) 
directly from these representations, (without having to calculate n*). 
This and other results will follow directly from some simple linear 
algebra. In Lemmas 1 and 2 below, we use the following observation: 
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if 3C* and 3” are two subspaces of an n-dimensional real linear space 3, 
such that 3C” and 3C” have no non-zero vectors in common, then 
dim (35¢*) + dim (#) S n. 


Lemma 1: Let D be any real n X n matrix (not necessarily symmetric), 
and let 3” be any subspace of 3C such that for all x ¥ 0, xe 3, we have 
(2, Dx) < 0. If W and Z are n X n matrices such that (x, Wr) = 0 
for allxe Hand D = W — Z, then rank Z = dim (8). 

Proof: Let 3* denote the null space of Z. Clearly (z, Dx) = 0 for 
ze 3c’, and dim (3¢*) = n-rank (Z). Also, 3¢* and 3C” have no non-zero 
vectors in common so that dim (3¢°) + dim (3¢*) = dim (3¢°) + 
n-rank (Z) S nor rank (Z) = dim (3). 

The following result is a dual version of the above: 


Lemma 2: Let D be any real n X n matrix and let K* be any subspace 
of 3 such that for all x # 0, xe 3C*, we have (x, Dr) > 0. If D= W - Z 
with (x, Zx) = 0 for all x, then rank W = dim (8 "*). 

The following corollaries are immediate from the two lemmas. 


Corollary: Let D be a real n X n symmetric matrix with n* and n™ 


positive and negative eigenvalues respectively, and le D = W — Z. If 
(x, Zz) = O[(x, Wx) = 0) for all x, then rank (W) = n° [rank (Z) = n}. 
Corollary: Let D be the distance matrix of a graph G. Then N(G) 2 
max {n*,n}. 


Proof: From condition 3 of Theorem 1, we have 2D = CCt — EEt 
where C and £ are n X N matrices. Then application of the previous 
corollary yields 


N 


IV 


rank C = rank CCt 


IV 
= 
+ 


and 
N 


IV 


rank # = rank EEt =n. 


Note that the proof of this result does not depend on the integer 
nature of the matrices involved. 

We mentioned in the previous section that factorizations of the form 
D = PP* have special significance. We now demonstrate this. 
Theorem 5: Let D be any symmetric non-zero matrix satisfying D = PP* 
for some binary-valued n X L matrix P. Then n™ = land n = 
rank (D) — 1. 


Proof: Let e, and e, be vectors of length Z and n respectively, all of 
whose elements are 1. Then e,e, is a rank one, L X n matrix of I’s 
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and P* = e,e, — P*. Thus D = Pe,e, — PP*. But PP is non-negative 
definite and rank (Peze,) = 1. Hence, Lemma 2 implies n* < 1. How- 
ever, D;; = 0 for all7 andj and D ¥ 0 so that e'De, > 0. Thus n* = 1, 
and since n* + n = rank (D), we have n” = rank (D) — 1. 

Incidentally, the above shows that every row of P has the same 
number of 1’s. 

This Theorem is useful because it enables one to determine n* and n™ 
from rank D. Thus, if a distance matrix of a graph G satisfies D = PP*, 
then N(G) = max (1, rank D — 1) = rank D — 1. This result has 
immediate application to complete graphs, trees, cycles, and cubes 
considered in the previous section, since in each of those cases, we 
exhibited a factorization of the form PP*. We remarked previously 
that for any graph that admits an address matrix H having no a’s as 
entries, D = PP* for some P. We now derive an upper bound on the 
number of rows of H that do not contain d’s. 


Theorem 6: Let D = PQ‘ be symmetric, and suppose that r rows of P 
are the complements of the corresponding r rows of Q. Thenr Sn—n +1. 
Proof: D = PQt = P(Qt+ Pt — Pt) = P(Qt+ Pt) — PPt 

From Lemma 2, 


rank (Qt + Pt) = rank (P(Qt+ P*)) 2 n°. 


But Qt + P*t has n columns, r of which have ones in all positions. Hence 
rank (Qt + Pt) S$ n — r + 1, and the desired result follows. 


Vv. AN ALTERNATE ADDRESSING SCHEME 


In this section we discuss an addressing scheme that achieves minimum 
distance routing, and is very simple to construct for all graphs (in fact, 
D could be any non-negative integer-valued matrix; this could often 
arise whenever preference for routes is not dictated by just the path 
length). Let s be the diameter of G, a graph on n vertices. We consider 
a PQ* addressing of G with P and Q matrices each of order n X ns. 
The addressing that follows is a simplified version of the construction 
used in Theorem 4. The 7th row of P, P; , is zero everywhere except in 
positions (s(¢ — 1) + 1) to sz where it is one. The zth row of Q, Q; , is 
constructed as follows: for every j, there are exactly D;,; 1’s in positions 
(sj — 1) + 1) to sj with the rest of the entries of Q; equal to zero. 
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S 2s 38 nS 
1 J _ 
1---1 0---0 0--:-0 0---0 
eet 266, Oe 
P= 
0 0 0 0 0---0 1 1| 
8 Dis 2s D3 3s Din ns 
———_—_——{ i a aa ~ 
J J d, L 
0 Qo see Tse GO: Reseed ea gee See ee 
(xi O° 6 or OU veep oo3.0 (Mee ere 
Q= 
foe", 168 A ee 6S OC) Tek el oaet§ 0 hase 0 
—Y SSeS SS 
Di Dis D,3 


Obviously, the length of addresses in this scheme in ns. The P matrix 
defined above, rows of which are addresses to be prefixed onto messages, 
is the same for all graphs with diameter s. The Q matrix, rows of which 
are stored in the loops in the network, contains the information that 
identifies a particular graph. But the simplicity of the rows of P can be 
exploited to minimize the length of the message addresses in the following 
way: since the integer 7 (Sn) uniquely specifies row P; , the message 
address need only consist of the binary representation of 7, which, 
of course, requires at most log, n bits. The set of binary representations 
of the integers 1, --- m is the set of addresses to be prefixed to messages. 
In each loop, a device is placed which generates row P; from the binary 
representation of 7. The distance calculation, say in loop 7 is accomplished 
by forming the scalar product of the generated P; sequence with the 
stored Q; sequence. This calculation can be mechanized by an “and 
gate’ followed by a counter (as can be any of the PQ* decompositions 
discussed in Section III). 
This scheme has the following advantages: 


(t) The simplicity of constructing addresses by inspection of the D 
matrix is important in the case of large graphs with no special structure. 
Even if one were to find good minimum N(@) and L(G), we suspect 
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that the construction of a minimal length addressing will be very 
complicated. . . 
(27) Since present plans for length of message blocks envisage lengths of 
perhaps a few thousand bits, the length of the message address is 
an important parameter in any large loop system. The method of this 
section guarantees minimum length addresses (log, ») provided that 
the graph is not constrained to have a particular structure. 


One of the disadvantages of the present scheme is that it requires a 
large amount of storage in each loop. However, this can be remedied 
to a certain extent by coding the number of 1’s in each cell of length s 
of the rows of Q in binary form with log, s bits. Thus the storage re- 
quirement can be reduced fron ns bits to n log, s bits. The engineering 
aspects of this scheme and some of its modifications will be the subject 
of a forthcoming memorandum. 


VI. ANOTHER ALTERNATE ADDRESSING SCHEME 


We present here a coding scheme of the type PQ‘, which can address 
all graphs, using ternary logic symbols +1, —1, 0, and requiring ad- 
dresses of length 2(n — 1) for a graph on n vertices. The scheme uses 
the fact that for any graph there exists a numbering of vertices such 
that every vertex 7 is adjacent to some vertex 7 where 7 < 2. This, of 
course, implies for some pairs 7 and j that | D;,,— Din | < 1 for all m. 
We exhibit the construction inductively. Let D(r) be the r X r sub- 
matrix of D consisting of the first r rows and columns of D. 

Suppose D(r) = P(r)Q*(r) where P(r), Q(r) are (+1, 0)-valued 
matrices of order r X 2(r — 1). Then it is easily verified that 


Q(r) Q(r) 


De +1) = [PO 0 “ At 0 
P(r) 1 0 ot 1 


Here P,(r) is the sth row of P(r) where the (r + 1)st vertex is adjacent 
to the sth vertex, or | Doriitym — Dem | < 1 for all m. The symbol 0 
denotes the r X 1 matrix of zeros, A; , the jth component of the r X 1 
matrix A is +1, —1 or 0 according as 


Do+iyi — Ds; = +1, —1, or 0 respectively. 


Calling the obvious matrices P(r + 1) and Q(r + 1), we have D(r + 1) 
= P(r + 1)Q(r + 1) where P(r + 1), Q(r + 1) are (y + 1) X 2r matrices. 
Since D(2) = [} {] [? 4], the construction is complete. Because of the 
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structure induced on P and Q by this construction, one can code the 
rows P and Q making the length of the code smaller, at the expense of 
increasing the complexity of mechanization. The way this address is 
constructed simplifies changing the address when a new node is added 
to an existing coded graph. Notice also that the positivity of the ele- 
ments of D is unnecessary for this construction. 
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On the Accuracy of the Depletion Layer 
Approximation for Charge Coupled Devices 


By J. McKENNA and N. L.SCHRYER 
(Manuscript received February 28, 1972) 


In this paper we examine the accuracy of a variant of the depletion 
layer approximation as used in the analysis of the electrostatic potential, 
in the absence of mobile charge, in various charge coupled devices. The 
approximation is a linearization of the nonlinear potential equations 
and is wmportant in the numerical solution of two-dimensional problems. 
A one-dimensional model 1s solved both exactly and by means of the depletion 
layer approximation and the two results are compared. We conclude that 
this variant of the depletion layer approximation is excellent for buried 
channel devices, and adequate though not as good for surface devices. 
Some criteria are given which indicate the excellence of the approximation 
and can be used to estimale errors in two-dimensional calculations. 


The purpose of this paper is to discuss the accuracy of a variant of 
the well known depletion layer approximation’’” used in the analysis 
of the electrostatic potential, in the absence of mobile charge, in various 
charge coupled devices (CCD’s).*'* The approximation is a linearization 
of the nonlinear potential equations and is important in the numerical 
solution of two-dimensional problems. An application of these results 
is given in Refs. 5 and 6. We consider a one-dimensional model, which 
can be solved exactly, and examine in some detail the errors involved 
in the approximation. This analysis is a generalization of similar results 
in a previous paper.’ 

Consider a buried channel CCD’ formed by a first layer,O < y Sh, 
of SiO, ; a second layer, h; S y S hz , of completely ionized p-type 
Si; and a third layer, h. S y < ©, of uniformly doped n-type Si. (See 
Fig. 1) Then the dimensionless equations governing the potential 
g(x, y) are 


Vela, y) = 0, 0< y< h, ’ (la) 
V'ex(x, y) = oly), h<y<h, (1b) 
V’'e3(t, y) = exp (y.lz,y)) -—1, bh<y<-, (1c) 
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V2 baly)=a (y) 


V7b3 (y)=exP[3(y)]-1 n—TYPE Si 


DEPLETION REGION 





Fig. 1—A schematic diagram of a buried channel CCD. 


where ¢;(z, y) is the value of ¢(z, y) in the j-th layer and V’ = 0°/dz” + 
0° /dy”. In eq. (1), the dimensionless potential y(z, y), and distances, 
x, Y, hy , he , are related to the dimensional quantities ¢*(x*, y*), 2*, 
y*, h*, h*, respectively, (measured in MKS units) by 


et, v) = pA ete, y*), (2a) 


t= a*/Xp ’ y= y*/Xp ) hy = h¥/Xp ’ hy = h¥/Xp ) (2b) 
where the Debye length Xp is defined by 
An = (ekT/e’N¥)', (2c) 


and e > 0 is the electronic charge, k is Boltzmann’s constant, 7’ is 
the absolute temperature, e, the permittivity of the Si, and N¥ is the 
number density of the donors. Also, the dimensionless density of 
acceptor ions in the p-layer, c(y), is related to the physical density 
Ni(y*) by 


o(y) = Nik(y*)/N3§ . (2d) 
We assume that at the boundary y = 0, ¢,(z, 0) is completely specified. 
gilz, 0) = Vo(z) = 0, (3a) 


while at the base y = ©, 
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gs(z, ©) = 0. (3b) 


At the internal boundaries, y = h, , y = h2, ¢ satisfies standard electro- 
magnetic boundary conditions: 


0 0 
gi(x, hy) = 92(x, h), n Gi (v, hi) = an (z, hi) + Q,, ; (4a) 
£653) SuG LL iy eee ae (4b) 
Ze» 1b2 BX) Fl2)5 dy “°°? dy “0%? 


where 7 = ¢€,/€ is the ratio of the permittivities of SiO, to Si, and Q,, 
is a possible trapped charge at the oxide-semiconductor interface. The 
remaining boundary condition is taken care of by assuming periodicity 
in the x-direction 


e0,y) = oL,y), ODSy¥< 0 (4c) 


(see Ref. 5 for a fuller discussion of these boundary conditions). 

In the buried channel CCD, the parameters are adjusted so that 
the potential minimum occurs in layer 2, h; S y S hz , where the mobile 
charges are to be trapped. For most purposes, it is only necessary to 
know the potential accurately in layer 2. The doping in the p-layer is 
typically introduced by ion implantation followed by drive in diffusion.* 
Then the acceptor density can be given quite accurately by 


o(y) = C, exp ‘-(2 ah) In . —1 (5) 


where c, is the surface density of charge due to the original ion im- 
plantation. In some applications it 1s convenient to replace o(y) by 
its average 





Vai 6 


1 ha 
dy = 
he =) [ oly) 2 ~VIne, 


where erf is the error function.” In a typical buried channel CCD, 
region 1, the oxide layer, is 0.lu thick and y = 1/8. Region 2, the 
p-layer is 5y thick with an initially implanted charge density at the 
surface of 4.6 X 10'°/em’*, which corresponds to an average doping 
of ~2 X 10'°/em*. Region 3, the n-layer, has a doping of 10'*/cm’. 
In our dimensionless units, this corresponds to h,; = 0.24, h, = 12.24, 
c, = 46 and « = 19.71 at room temperature. (In all the examples 
where ¢ is used, we take ¢ = 20 for simplicity). The potential has been 
scaled so that 1 volt corresponds to 40 dimensionless units. Typical 
plate voltages are in the range —400 < V, S 0. Although we shall 








erf(WInc,) —1, (6) 


¢= 
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typically assume Q,, = O in this paper, we include it because it is 
necessary in some applications. 

If layer 2 is removed, we have a surface CCD.’ In a surface CCD, 
the mobile charge is trapped at the oxide-semiconductor interface, and 
only near this interface does the potential need to be known accurately. 
The remaining dimensions are typically the same as for a buried channel 
CCD, and the plate voltage is in the range —400 $ V, S —160. 

If 3(z, hz) is large and negative for 0 S x S JL, then for y > A, 
and near hz , exp (¢;(2, y)) « 1 and the equation for ¢; is approximately 
V'¢3; = —1. The region where this holds is called the depletion region, 
and at any point z, the depth of the boundary of this region below 
y = hz , R(z), is called the depletion depth (see Fig. 1). For (a, y) such 
that y > R(x) + he, | ¢s(x, y) | <« 1 and we have approximately Vg. = 
¢2 . In most problems R(x) is not a constant, but it is nevertheless 
often a good approximation to replace R(x) by some suitable constant 
R. This suggests replacing the system of eqs. (1) by the system of 
linear equations 


V'*vilz, y) = 0, O0<y<h, (7a) 
VW" pola, y) = oy), h<y<h, (7b) 
V'vs(a, y) = —1, he<y<h=h+R, (7c) 
V'vle, 9) = vAlt, 9), n= in + R<y<%, (7d) 


(See Fig. 2). In addition to ¥, , ¥. and ys; satisfying boundary conditions 
(3a), (4a), (4b) and (4c), we have the boundary conditions 


W(x, hs) a Walz, hs), a (x, hs) = oa ; hs), W(x, oo ) = 0 (8) 
which must hold for 0 S$ « S L. 

This paper ts devoted to studying how the choice of the pseudo-depletion 
depth R affects the accuracy with which W(x, y) approximates ¢(z, y). 
We do this by considering the special case Vo(x) = Vy < 0, where Vo 
is constant, in boundary condition (8a). For this choice of Vo(z), eqs. 
(1) and (7) reduce to systems of ordinary differential equations the 
solutions of which are functions of y only. In all that follows, we drop 
all reference to x, replace V” by d’/dy” in eqs. (1) and (7), and consider 
only the one-dimensional problem. The remainder of the paper is 
devoted to calculating and comparing the values of ¢.(h2) and g$(he) 
with ¥.(h.) and 5(h.), where’ = d/dy. It is easy to show that for 
hsysh, 
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V? baly)=a(y) p-TYPE Si 


Vewaly)=-1 


DEPLETION REGION 


n—TYPE Si 





Fig. 2—A schematic diagram of the depletion layer approximation to a one- 
dimensional, buried channel CCD. 


ley) — ¥y) | Ss | ga(h2) — Wolhe) ; 
le’(y) — WY) 1 S | ¢2(he) — Pathe) |- (9) 


We only give the analysis for the buried channel CCD, that is very 
similar for the surface CCD. We give numerical results for both. 
We begin by solving for g.(hz) and ¢}(hz). We can now write 


gily) = ay + Vo, (10) 
ese) = [ & ~ Wo a + Wy — hy) +e (11) 

In addition, a first integral of (1c) subject to (8b) is’ 
3{es(y)}* = exp {ya(y)} — gay) — 1. (12) 
Making use of the boundary conditions (4), it is now easy to show that 
a=1(o-[' wa+e.), (18) 


c= (m—-i + bem Ho.- [ (e-m +2 )o@ a, a9 
u] 0 ha 0 


ip? =e —e—1. (15) 
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Equations (14) and (15), for c = g(hz) and b = (hz), can be solved 
numerically. In fact, for most cases of interest c «< —1, so the term 
e° in (15) can be neglected, yielding 


b = gi(he) = —(ha -h+ ts) 


tyl(i—n +4) ~2] vi41- if (nt ew aco,,| , (16) 


n ” 
C= (hz) = —ib° — 1. (17) 
It can be shown that for a surface CCD, v3(h,) and ¢2(h1) can be ob- 


tained from (16) and (17) by setting o(y) = O and hi = h.. 
To solve eqs. (7), we can write 


wily) = Ay+ Vo, Osysh, (8a) 
H0) = [@- Wo) +BY) FC, bSy Sh, (Bb) 
vay) = —2y — he)’ + Dy — he) +E, ho Sy Shs, (18) 
valy) = de" Rp ah t+ Rsy < @, (18d) 
where the boundary conditions at y = 0 and y = @ have already 


been used. Upon employing the remaining boundary conditions at 
y = h,, hg and hy + R = A; , we obtain six equations relating the seven 
constants A, B, C, D, E, 6 and R. 

Our main interest is in cases where R is given, so that A, B, C, D, E 
and 6 can be determined in terms of R. However, we first consider the 
case where 6 is specified. This will provide useful information about 
picking an optimal value of R later on, and will also give a useful 
criterion for estimating how good the approximate solution is when R 
is specified. 

In Ref. (1) the effects of requiring that 6 = 0 (the usual depletion 
layer approximation) or 6 = —1 were studied, and it was shown that 
the choice 6 = —1 yields a better approximate solution. It should be 
noted that requiring 6 = —1 is equivalent to requiring that Y3’(hs) = 

4/(hs). We show here that if the value of 6 is specified, the choice 
6 = —1 does yield the best approximation, but the values of p.(he) 
and ¥i(hz) are relatively insensitive to the choice of 6. It is easy to 
show that in terms of 6, 
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A hy, h, . 2 
ae eee em ae ae See + 6 + 26 — 2V, 


+ of ( digas )o( dé — 2m Q, ay (19a) 


a i 
A=—R +4(a. — [ a(t) dé — s) (19b) 
Wh) =B=D=R-—6, (19¢) 
Y(los) =C =E = 64+ Rs — 4R?. (19d) 


For a surface CCD, the expressions for ¥3(h,) and y.(f,) are given 
by (19) when o(y) = Oand hy = h.. 

As a first example, pick the typical parameters mentioned earlier, 
h, = 0.24, hp = 12.24,Q,, = 0,7 = 1/3, Vo = 0, and let o(y) = « = 20. 
In Fig. 3, using eqs. (16), (17) and (19), we show a plot of 

e(y) = | Wolhe2) — ¢e(he) [/| e2(he) | (20a) 


and 
e(y’) = | W2(h2) — oo(he) I/| ga(he) | (20b) 


as functions of 6 for —10 S 6 S 10. It can be shown from (17) that 
these relative errors have a minimum at 6 = —1 (and are even functions 


e(w’) 


e(ys) OR e(p’) 


1073 


12 —8 —4 0 4 8 12 
5 


Fig. 8—Plots of e(y) and ev’) for a buried channel device with h,; = 0.24, 
he = 12.24,¢ = 20,Q;; = 0,7 = 3, and Vo = Oas functions of for -—10 S$ 6 S 10. 
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of 6 + 1), but over the whole range of 6 shown, the maximum relative 
error in W2(h2) is only 1.27 percent and in W4(h2) is only 2.28 percent. 

As a second example, consider a buried channel CCD with the same 
parameters as in example one, except that in this case we assume that 
a(y) is given by eq. (5) with c, = 46. In Fig. 4, we show a plot of e(y) 
for this case as a function of 6 for —10 S$ 6 S 10. The results are 
qualitatively the same as in the first example, except that the relative 
error is about twice as large. It should be remarked that the potential 
distributions in the first and second examples are quite different. The 
potential minimum in example one occurs much nearer the bottom of 
the p-layer than in example two. In example one, g(h2) = —1034.5 
while in example two, ¢(h2) = — 530.6, so we would expect the depletion 
layer approximation to be better in example one than in example two, 
as it is. 

As a third example, consider a surface CCD with h, = h. = 0.24, 

= 0,Q,, = 0, 7 = 1/3, and V, = —160. In Figs. 5a and 5b, using 
eqs. (16), (17) and (19), we again show plots of e(y) and e(y’) as func- 
tions of 6 over —10 < 6 < 10. We see a maximum relative error of 
1.58 percent in y.(h,) and 18.3 percent in W3(h,). 

Since the approximation is based on the assumption that ¢(he) is 
large and negative, the magnitude of y(h.) must affect the accuracy 
of the approximation. For the buried channel CCD, g(hz) is S$ —500 
for all Vo S 0, so in this case, the approximation is very accurate for 


a ie 


10-4 
-12 -10 -8 -6 -4 -2  O 2 4 6 8 10 12 


Fig. 4—Plot of e(y) for a buried channel device with hy = 0.24, hp = 12.24, 


cs = 46,Q;, = 0, » = 4, and Vo = 0 as a function of 6 for —10 S$ 6 S 10. 
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Fig. 5(a) and 5(b)—Plots of e(y) and e(y’), respectively, for a surface device with 
= h, = 0.24, 6 = 0, Qss = 0, 9 = 3 and Vo = —160 as functions of o for 
~10 $6 < 10. 


any negative V, . However, in the surface CCD, as V, — 0 we have 
e(h,) — 0. Thus for “small” V, , we expect the approximation to be 
quite bad. Figure 6 gives a plot of the relative error e(y) for a surface 
CCD as a function of V, using the same parameters as given for Fig. 5 
and 6 = —10. To maintain an error of the order of 1 percent for all 
6 = [—10, 10], it is clearly necessary to keep Vp S —160. 

We turn now to the main problem where R is arbitrarily specified. 
This occurs typically in the numerical solution of two-dimensional 
problems where it can be reasonably assumed that the depletion depth 
R(x) does not vary too much with the various plate potentials. Then 
an average depletion depth is estimated, sav from eq. (19a) by setting 
6 = —1 and choosing for V, some average plate potential. 

We now investigate the errors involved in this approximation. If 
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0 —160 —320 —480 —640 —800 —960 
(0) (=Vo) 


Fig. 6—A plot of e(y) for a surface device with hi = hz = 0.24, 0 = 0, Q;; = 0, 
nm = 4and 6 = —10 asa function of Vo for —800 S$ Vo S —40. 


R is given, then it is easy to show that 


vite) =B=d=|["(e-n+™)@a-Ba.-v 


+H +k] /[Rtitn—n +4], (21a) 

_1 . ) 
ae (3 — [ af) dé + Que) , (21b) 
Wille) = C = E = $R?+ R—- (R+ DB, (21c) 
5=R-B (21d) 


As before, the above expressions give y¥2(h,) and y3(h,) for a surface 
device when we set o = Oandh, = h,. 

As a fourth example, we consider a buried channel CCD with h, = 
0.24, ho = 12.24, o(y) = @ = 20,Q,, = 0, 7 = 3, and we choose R = 
50.98. From (19a), this corresponds to the “‘true’”’ depletion depth R 
(with 6 = —1) for a plate voltage V, = —400. In Fig. 7, using eqs. 
(16), (17) and (21), we give plots of e(W) and e(y’) as functions of Vo 
for —800 < V, S 0. In Fig. 8, we give a plot of 6 as a function of Vo 
from (21d) for the same data, and a plot of the “true” value of R 
from (19a) for 6 = —1 as a function of V, . Several important con- 
clusions can be drawn from these graphs. First, over —800 S$ V, S 0, 
¥(y) is an excellent approximation to ¢(y) for h, S y S h, , the maxi- 
mum relative errors in Y2(h2) and W5(h2), which occur at V, = 0, being 
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10-2 
> 
rot e(y') 
fom -3 
e10 
= 
® e(yp) 

1074 

0 —160 —320 —480 —640 —800 -960 


$ (0) (=Vo) 


Fig. 7—Plots of e(y) and e(y’) for a buried channel device with h; = 0.24, 


ue = 12.24,¢ = 20,Q;, = 0,7 = 4 and R = 50.98 as functions of Vo for —800 S 
os 0. 


0.41 percent and 0.74 percent respectively. Second, y(y) approximates 
g(y) excellently over h, S y S h,., even though a relatively large error 
has been made in estimating the true value of R. In fact, R = 50.98 


differs from the ‘‘true” value by 11.6 percent at Vo = —800 and by 

8 58 

6 56 
4P 54 

2 52 

e 0 50 R 

-2 48 
~4 46 
-6 44 
-8 42 
40 

0 - 160 -—320 -— 480 -640 - 800 — 960 

(0) (=Vo) 


Fig. 8—Plots of 6 and “true” & for a buried channel device with h, = 0.24, 
- = 12.24,¢ = 20,Qs, = 0,7 = 3 and R = 50.98 as functions of Vo for —800 S 
<0 


o= . 
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0 -160 -320 —480 —640 =800 —960 
$ (0) (=Vo) 


Fig. 9—Plot of e(y) for a buried channel device with hi = 0.24, he 


12.24 
Ce = 46,0;; = 0,7 = 


4+ and R = 39.66 as a function of Vo for —800 <= Vo = 0. 


12.7 percent at Vo = 0. Finally, an estimate of the error in ¥(h2) can 
be obtained from monitoring 6. For Vo S 0, as long as —7 S 6 S 5, 


the error in ¥2(h.) and (hz) is less than 1 percent. 
Next, consider a buried channel CCD with h, = 0.24, h, = 12.24, 


46 
4 42 
° 0 38 i 
-4 34 
-8 30 
0 — 160 —320 -—480 —640 —800 — 960 


(0) (=Vo) 


Fig. 10—Plots of 6 and ‘‘true” # for a buried channel device with hi = 0.24, 


he = 12.24, c,; = 46,Q,; = 0,7 = 3 and R = 39.66 as functions of Vo for —800 $ 
Vo S 0. 
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Fig. 11(a) and 11(b)—Plots of e(y) and e(y’), respectively, for a surface device 


with hy = ho = 0.24, ¢ = 0,Q;; = 0, 7 = 4 and R = 26.55 as functions of Vo for 
—800 S$ Vo S —40. 


Q,, = 0, 7 = 3, with o(y) given by (5) with c, = 46, and with R = 
39.66. From (19a), this corresponds to a “true” value of R (with 6 = —1) 
for a plate voltage of V) = —400. In Fig. 9, we plot e() as a function 
of V, for —800 S V, S 0. In Fig. 10, we give a plot of 6 as a function 
of V, from (21d) for the same data, and a plot of the ‘‘true” value of R 
from (19a) for 6 = —1 asa function of V, . The results are qualitatively 
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the same as for example four, but as we should expect, the errors in 
using the depletion layer approximation are greater in example five 
than in example four. 

As a final example, we consider a surface CCD with hi = h, = 0.24, 
o = 0,Q,, = 0, = 4 and R = 26.55. From (19a), this corresponds 
to a “true” value of R (with 6 = —1) for a plate voltage of Vp = —400. 
In Figs. lla and 11b, using eqs. (16), (17) and (21), we give plots of 
e(y) and e(y’) as functions of V, for —800 $ V, S —40. In Fig. 12, 
we give a plot of 6 as a function of V, from (21d) for the same data, 
and a plot of the ‘“‘true” value of R from (19a) (for 6 = —1) as a func- 
tion of V, . Hence, we note that over the range —800 < V, S$ —160, 
the error in y.(h,) remains remarkably small, less than 1 percent. 
However, ¥3(h;) does not approximate g3(h,) nearly as well, especially 
as V) — 0. Note also that in this case much larger errors have been 
made in estimating R, the error being 39 percent at V) = —160 and 
44,2 percent at V, = —800. Again, 6 provides a check on the accuracy 
of the approximation, with —10 S 6 S 8 assuring an error of less than 
1 percent in approximating ¢(y) inh; S y S h.. 

We can conclude that this variant of depletion layer approximation 
works excellently for the buried channel devices in providing approxi- 
mations to ¢(y) and ¢’(y) in the region h, S y S h, . In addition, the 
parameter 6 is a good index of the excellence of the approximation. 
If in a two-dimensional calculation this variant of the depletion layer 
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Fig. 12—Plots of 6 and “true” # for a surface device with hy = hz = 0.24, 0 = 0, 
Q:s = 0,7 = 4 and R = 26.55 as functions of Vo for —800 S$ Vo S —40. 
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approximation is used, then a value of the pseudo-depletion depth R 
is estimated, say from using an average plate potential in eq. (19a). 
Then, so long as the solution of this approximate problem remains 
in the interval [— 10, 10] along the line y = R + hz, one can feel reason- 
ably confident that the approximate solution differs from the true 
solution, in the region of interest, by at most a few percent. It should 
be noted that the depletion layer approximation is poor near y = 
ho + R, but in many important cases this is unimportant. Finally, 
for surface CCDs, the approximation is good to g(y) but not as good 


to ¢’(y). 
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Perceptual Evaluation of the Effects of 
Dither on Low Bit Rate PCM Systems 


By L. R. RABINER and J. A. JOHNSON 
(Manuscript received April 10, 1972) 


It has previously been shown that by adding a pseudo-random ‘‘dither”’ 
noise to a signal to be quantized, and by subtracting an identical noise 
sequence from the quantizer output, it is possible to break up undesirable 
signal-dependent patterns in the quantization error sequence without 
increasing the variance of the error. The effect of the dither noise becomes 
significant when the number of bits per sample is less than about six. An 
experimental evaluation of the perceptual effects of dither on speech has 
shown: 


(t) strong preferences for dithered speech over straight PCM encoding 
at identical bit rates, 

(iz) for low bit rates (2-4 bits/sample), a preference for dithered speech 
over PCM encoded speech even when the PCM speech had one more 
bit per sample than the dithered speech, 

(122) an increase in word intelligibility for dithered speech over PCM 
speech when 4 to 6 bits/sample were used, 

(wv) a decrease in word intelligibility for dithered speech over PCM 
speech when 2 to 3 bits/sample were used. 


I. INTRODUCTION 


When a signal, such as a speech waveform, is quantized, the quantiza- 
tion error waveform is usually correlated with the original signal. This 
correlation is virtually inperceptible when the quantization is quite 
fine-l.e., a large number of bits/sample. For crude quantizations, 
however, the correlation becomes quite large and the quantization error 
is easily perceived. As a result, it can become quite disturbing to listen 
to speech quantized to a low number of bits/sample for an extended 
period of time. In such cases, techniques that decorrelate the quantiza- 
tion error from the signal are attractive, even if they do not increase 
the signal-to-noise ratio of the system. Dithering is such a technique in 


1487 
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which a pseudo-random “dither” noise is added to the speech before 
quantizing, and then the identical noise is subtracted producing a 
quantization error which is uncorrelated with the original speech wave- 
form.’ Figure 1 shows a comparison between straight PCM and a 
system in which dithering is used. In an earlier work, Jayant and 
Rabiner’ discussed several theoretical issues involved with dithering 
and demonstrated its utility for the quantization of speech signals. In 
this paper, we present experimental results on the perceptual effects of 
dither on both the preference and intelligibility of PCAI encoded speech. 


II], PREFERENCE EVALUATION TEST 


The purpose of this experiment was to determine the perceptibility 
of the decrease in correlation between the quantization error and the 
original speech, as a function of the number of signal bits. 

The stimuli used in the preference test were a set of ten sentences 
chosen from a list of ‘everyday speech” sentences® compiled at the 
Central Institute for the Deaf. The sentences used are shown in Table I. 
These ten sentences were spoken by a General American speaker, 
digitized at a 10 kHz rate with 16 bits/sample, and stored on the disc 
of the DDP-516 computer. 

In order to limit the number of stimuli to be used in the paired- 
comparisons preference test, the number of bits/sample was restricted 
to the range of 2 to 6 bits. Therefore, there were ten distinct stimuli in 
the test, i.e., (five possible values for the number of bits) X (two types 
of quantization—dither or straight PCM). For notational convenience, 
the stimuli were coded using a two-digit code. The first digit refers to 
the number of bits/sample (i.e., 2-6) and the second digit specifies the 
type of quantization. A O in the second digit means straight PCM 
encoding, whereas a 1 in the second digit means dithered speech. Thus 
stimulus 31 has 3 bits/sample and uses the dither noise, whereas stimulus 
50 has 5 bits/sample and does not use dither noise. 

Since there were ten distinct conditions to be evaluated, a complete 





N 
X L—BIT Ma wea 
QUANTIZER QUANTIZER | (X+N)o 
STRAIGHT PCM DITHER SYSTEM 


Fig. 1—Block diagrams of a straight PCM system and a dither system. 
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TaBLe I—SentTENCES USED IN PREFERENCE TEST 








1. Walking’s my favorite exercise. 
2. Here’s a nice quiet place to rest. 
3. Our janitor sweeps the floor every night. 
4, It would be much easier if everyone would help. 
5. Good morning. 
6. Open your windows before you go to bed. 
7. Do you think she should stay out so late. 
8. How do you feel about changing the time when we begin work. 
9. Here we go! 
10. Move out of the way. 


paired-comparison preference test involved 100 pairs. These 100 pairs 
were randomly generated by a DDP-516 program which randomly 
accessed each of the ten stimulus sentences ten times in the course of 
the experiment. Each of the 100 stimulus pairs was recorded on magnetic 
tape for offline running of the experiment. 

Ten subjects participated in the experiment. Each subject was given 
the following instructions: 


“In this test you will be listening to pairs of sentences. Each of the 
two sentences (first is called A, second B) was processed by some 
type of speech transmission system. After you hear both sentences, 
there is a five-second interval in which you are to write down the 
sentence, A or B, you prefer, i.e., the type of transmission system 
you would prefer listening to for an extended period of time. You 
must choose either A or B-even if you have no preference.”’ 


The preference test required two 15-minute listening sessions per subject 
and was run on tio separate days. 


III. RESULTS OF PREFERENCE TEST 


For each of the ten subjects, a matrix of preferences was determined 
in which a 1 in a particular cell of the matrix denoted that stimulus B 
is preferred to stimulus A, and a 0 indicated the reverse condition. 
Table II shows the matrix obtained by summing the matrices for the 
ten subjects. Careful inspection of this matrix shows a strong preference 
for dithered speech over straight PCM encoding at a fixed number of 
bits/sample, and, in many cases, a preference for dithered speech at 
L bits/sample (L = 2-4) over straight PCM encoded speech at (L + 1) 
bits/sample. 

To verify these preference results, the data was analyzed using a 
multidimensional preference program of Carroll.* The program indicated 


1490 THE BELL SYSTEM TECHNICAL JOURNAL, SEPTEMBER 1972 


Taste [JI—Matrix or SuM OF PREFERENCES FOR PAIRED 
COMPRESSION PREFERENCE TEST 
Stimulus B 


20; 3! 8] 8 | 10} 10] 10 | 10} 10} 10] 10 





NS) 21 0 5 5 9 7/10;10/] 10) 10 9 
i 361 con sel oo: orl 40" doses a0 | 07 
: aril vor od a | Ol ao; ab ao aol ao: 
: io. | 2 aa es eo oe Oo tae 
: idl tole) oko al oa cet al aos 
A 50| o| of 2! 4/ 1/1 61 41/10/10) 10 


51]; O; O} OF 1 1/ 4] 0O|] 7] 8] 7 


60; 0; O}; OF} OF O 1 0!) 7 2) 3 
61 


oO 
fo) 
o 
o 
o 
wo 
me 
aS 


Number of preferences of B over A in 10 trials 


that the preferences were essentially one-dimensional (over 95 percent 
of the variance was accounted for by one dimension), and produced a 
graphical interpretation of the overall preferences which is shown in 
Fig. 2. Since the preference judgments were one-dimensional, all the 
conditions lie on a line. The direction of preference goes from left to 
right in terms of decreasing preference. Figure 2 clearly shows: 


(c) For a fixed number of bits/sample the dithered speech samples 
are always preferred to straight PCM encoding, 

(11) For 2-4 bits/sample, dithered speech is preferred to straight 
PCM encodings even with one extra bit/sample, i.e., condition 41 
is preferred to condition 50, condition 31 is preferred to condi- 
tion 40, and condition 21 is preferred to condition 30. 


Thus in some perceptual sense, dithered PCM speech has a one-bit 
advantage over straight PCM encoding under certain conditions. This, 


* 6160 51 41 50 31 40 2130 20 
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PREFERENCE DECREASING PREFERENCE 


Fig. 2—Ordering of the stimuli in terms of preference. 
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of course, is not correct in terms of physical measures such as signal-to- 
noise ratio, or, as we will see, word intelligibility. 

A complete analysis of variance was performed on the preference data 
and the results of this analysis are shown in Table III. The three factors 
and the number of levels of each are: 


(7) number of bits/sample (5) 
(tz) type of quantization (2) 
(177) subjects (10) 


The analysis reconfirms the conclusions already discussed in that the 
most significant effects (significance > 0.999 level) were number of 
bits/sample, and type of quantization. Subjects were significant at the 
0.95 level, and the interaction between bits and dither was also signi- 
ficant at this level. 


IV. WORD INTELLIGIBILITY TEST 


The purpose of the intelligibility test was to determine the effects 
of dithering on the intelligibility of isolated monosyllables. As discussed 
earlier, the effect of dither is to make the quantization noise act like 
an additive wideband uncorrelated noise. Earlier studies’ have indicated 
that such a noise tends to mask consonants, thereby lowering intel- 
ligibility. The effect of the correlated quantization noise on straight 
PCM encoding on word intelligibility was also measured. 

In this experiment, 200 PB words’ (Lists 2, 4, 5 and 6 in Ref. 5) were 
recorded, digitized, and stored on the disc of the DDP-516. The words 
were accessed at random, in groups of 50 (i.e., an entire list was processed 
before a new list was used), by one of the ten systems used in the 


TaBLE JIJI—ANALYSIS OF VARIANCE OF PREFERENCE DATA 














Degrees of Mean Significance 
Factor Freedom Square F-ratio Level 

Subjects (s) 9 1.8 2.3 0.95 
Type of 

quantization (TQ) 1 53.3 67.7 >0.999 
Number bits per 4 142.9 181.5 >0.999 
sample (NB) 

S x TQ 9 0.8 1.1 N.S.* 
NB X58 36 1.0 1.2 N.S. 
NB X TQ 4 2.3 2.9 0.95 
Residual 36 0.8 





* N.S. = not significant above 0.90 level. 
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preference test. The 200 words were divided into two tests of 100 words, 
each test containing 10 versions of each stimulus condition. The same 
ten subjects were used in the intelligibility test as in the preference test. 
The two tests were given on separate days to all ten subjects. 


V. RESULTS OF INTELLIGIBILITY TEST 


Table IV shows the average error scores as a function of the number 
of bits/sample, and the type of quantization. (The notation of the 
previous section is used again here.) These data are averaged over 
subjects and tests. This table shows that at 2 bits/sample, the PCM 
system has an error rate of 59.5 percent as opposed to 76 percent for 
the dither system, i.e., a decrease of 16.5 percent in word intelligibility 
due to consonant masking. At 3 bits/sample, the PCM system has an 
error rate of 34.5 percent whereas the dither system has an error rate 
of 46.5 percent. Thus even at 3 bits/sample, the masking of the dither 
noise reduces word intelligibility by about 12 percent. At 4-6 bits/sample, 
the dither system has lower error rates than the PCM system-the 
differences being 10 percent at 6 bits/sample, 1.5 percent at 5 bits/sample 
and 0.5 percent at 4 bits/sample. Thus only at 6 bits/sample is the 
error rate difference significant. The data of Table IV are plotted in 
Tig. 3 to show how the error rate varies with the number of bits/sample 
for the two systems. 

A complete analysis of variance was performed on the raw data of 
the intelligibility test. The four factors used in the analysis (and the 
number of levels of each factor) were 


(c) number of bits/sample (5) 
(27) type of quantization (2) 
(227) subjects (10) 
(iv) repetitions (2) 


TABLE IV—Worp Error Scores AVERAGED OVER 
SUBJECTS AND REPETITIONS 

















Error Rate 
Number of Bits |——————__——_—- 

per Sample PCM Dither Difference 
2 59.5% 76% —16.5% 
3 34.5% 46.5% —12% 
4 29.5% 29% 0.5% 
5 25% 235% 1.5% 
6 16.5% 6.5% 10% 
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Fig. 3—The percentage error for word intelligibility as a function of the number of 
bits/sample for straight PCM and dither systems. 


The results of the analysis are shown in Table V. The most significant 
factor was, of course, the number of bits/sample. The next most signifi- 
cant factors were subjects, repetitions, bits/sample X type of quantiza- 
tion, and bits/sample X repetitions. These results indicate a fairly 
large amount of learning between repetitions 1 and 2, as well as a lack 
of consistency between the intelligibility scores of the different subjects. 
VI. CONCLUSIONS 

The results of the preference test were quite encouraging in that 
subjects uniformly showed strong preferences for dithered speech over 
straight PCM encoding at all bit rates employed in the experiment. 
At the lower bit rates, the preference for dithered speech over higher 
bit rate PCM encoded speech presents strong evidence for the per- 
ceptibility and annoyance of highly correlated quantization noise. 

The word intelligibility tests showed that the wideband uncorrelated 
dither noise tended to mask the consonants more than the correlated 
PCM noise thereby reducing word intelligibility by about 14 percent 
at low bit rates. At the higher bit rates used in the experiment, there 
was no decrease in word intelligibility for the dither system, and, in 
fact, at 6 bits/sample, the dithered words were 10 percent more intel- 
ligible than the straight PCM encoded words. Since the average per- 
centage correct for the PCM system was 83.5 percent, an increase of 
10 percent is a significant increase in intelligibility. 

Overall, these experiments indicate that the use of dither noise in 
the range of 4-6 bits per sample has many beneficial effects. 
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TABLE V—ANALYSIS OF VARIANCE OF INTELLIGIBILITY DATA 














Degrees of Significance 
Factor Freedom Mean Square F-ratio Level 

Repetitions (R) 1 29.6 30.5 >0.999 
Subjects (S) 9 5.7 5.8 >0.999 
Type of quanti- 

zation (TQ) 1 5.4 5.6 0.975 
Number bits per 

sample (NB) 4 180.3 185.6 >0.99 
RxXS 9 2.1 2.2 0.95 
TQ xR 1 0.05 0.05 N.S.* 
TQ xS 9 1.9 1.9 0.90 
NB XR 4 TA 7.3 >0.999 
NB XS 36 Led 1.7 0.90 
NB X TQ 4 11.6 11.9 >0.999 
TQ xRXS 9 0.5 0.5 N.S. 
NB XRXS 36 1.1 1.1 N.S. 
NB Xx TQ xX R 4 0.3 0.3 N.S. 
NB x TQ XS 36 0.8 0.8 N.S. 
Residual 36 1.0 








* N.S. = not significant above 0.90 level. 
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Methods for Designing Differential 
Quantizers Based on Subjective 
Evaluations of Edge Busyness 


By J.C. CANDY and R. H. BOSWORTH 
(Manuscript received February 18, 1972) 


This ts a study of the visibility of television noise and its dependence 
on the instantaneous rate of change of the video signal. Novtse added to a 
picture tends to be least noticeable in regions where the brightness changes 
rapidly, but the relationship between visibility-of-noise and slope-of-the 
signal 1s dependent on the scene being displayed. Measurements are 
presented for four different scenes, and these data are used to design com- 
panding laws that minimize the visibility of noise from Differential 
Quantizers. The designs agree with those that have been satisfactory in 
practical applications. 


I. INTRODUCTION 


Experience with Differential Quantizers’'’’* teaches us that noise 
superimposed on brightness boundaries of a picture is less objectionable 
than similar noise on areas of uniform brightness. This work is an 
attempt to evaluate the phenomenon and make efficient use of it when 
designing quantizing scales for Differential Coders. We know that there 
is advantage in companding the quantization levels so that the largest 
steps are used for reconstructing rapidly changing signals and the 
smallest steps for slowly changing signals. In the past, the companding 
laws for differential quantizers have usually been obtained experi- 
mentally by trial and error, but these techniques are inadequate when 
there are a large number of levels. A more objective method is described 
here; it provides a basis for theoretical synthesis and evaluation of 
companding laws independent of the number of levels that are required. 

The essence of the method is a subjective evaluation of the dependence 
of the visibility of noise on the rate of change of the signal. The measure- 
ment is used to determine how the visibility of the net quantization 
noise in a picture depends on the companding law. Laws that minimize 
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this visibility are then obtained and they agree very well with those 
that have been obtained experimentally. Reference 4 describes a design 
of differential quantizers that is based on the probability of slopes 
occuring in pictures and on a weighting function for the visibility of 
noise. 


Il. EVALUATING EDGE BUSYNESS 


2.1 Simulating Edge Busyness in Pictures 


Noise added to a video signal at times when it changes rapidly has 
the appearance of “‘busyness”’ near brightness boundaries of the scene.’”® 
In order to measure its visibility, the noise must be introduced into 
the signal in a controlled fashion. A convenient method compares the 
rate of change of the signal amplitude with a threshold value, and adds 
calibrated noise at times when the threshold is exceeded. 

Figure 1 illustrates the circuit used for the experiment. The input, 
a signal that simulates that of a Picturephone® system, passes through 
two parallel paths. In the upper path the signal is first differentiated, 
and then full-wave rectified to obtain a voltage that is proportional 
to the magnitude of the rate of change of the input. Whenever this 
voltage exceeds a threshold value, a trigger circuit fires and gate G 
conducts. A random signal is then fed through a high-pass filter and 
added, under control of switch S, to the video signal. In the picture, this 
noise resembles the edge busyness which is characteristic of differential 
coding. Its amplitude is controlled by calibrated attenuator Al. 

The observer has the choice of two pictures by means of switch S: 
in position A, edge busyness is added to the scene; in position B, white 
noise is added to the entire scene. The observer varies the setting of 
attenuator B1 until, in his judgment, the quality of the picture is the 
same for the two positions of the switch. The experiment is repeated 
for various settings of the threshold and of the attenuator Al. 

Attenuators Al and Bl were calibrated to read the ratio of the 
rms amplitude of the noise to the peak signal. The rms value of the 
“edge-noise’’ was measured when the threshold was set at zero, in 
order that the noise be present at all times. This gives a measure of 
the noise in areas where it is displayed as opposed to a measure of the 
noise in the entire picture. 

The video signal used in the experiment was derived from a picture 
scanned with 271 interlaced lines at 30 frames a second. The visible 
portion of the picture was about 13 cm high and 14 cm wide. The peak 
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Fig. 1—The circuit used for simulating edge busyness and comparing it with 
white noise. 


luminance was about 70 foot lamberts and the room illuminance about 
60 foot candles. The white noise used in the experiments was approxi- 
mately Gaussian with a flat spectrum from 100 Hz to 1 MHz. The 
edge noise was obtained from a generator of random binary-words, 
driven at 2 MHz. The low-pass filters placed at the input and output 
of the system had frequency characteristics approximating a fourth- 
order Butterworth filter. At 1 MHz, the gains of the input and output 
filters had fallen by 6 dB and 15 dB, respectively. 

White noise with Gaussian amplitude distribution is used as a standard 
of comparison because techniques already exist” for describing the 
impairments introduced by uniform quantizing as an equivalent white 
noise. Moreover, it is a convenient and readily available reference 
for testing quantizers. A quasi-random binary signal is used as the 
source of edge noise because it is easily gated and may be synchronized 
to the scanning waveforms. 


2.2 Comparing Edge Busyness with White Noise 


The success of the experiments depends on the ability of observers 
to compare different impairments in a picture. These, and previous 
experiments’ have demonstrated that experienced observers have 
little difficulty in making the measurements. Inexperienced observers 
either tire or go through a learning period during which their measure- 
ments are not always reproducible. The results reported here were 
obtained by experienced observers. Their measurements tend to be 
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consistent and reproducible provided the picture quality is reasonably 
good; that is, having signal-to-noise ratios in excess of 40 dB. Eight 
observers took part in this study, and the differences between their 
opinions were small. The difference could usually be expressed as a 
variation of less than +1 dB in the sensitivity to noise. Most of the 
results presented are those obtained by one of the observers. His mea- 
surements have been compared with the other observers and are repre- 
sentative of this population. In future studies it may be useful to 
examine the differences between observers more carefully, but here 
we shall investigate the much larger differences that occur when scenes, 
and the type of noise is changed. 

Three main experiments are described, their purpose is to find the 
following: 


(t) The relationship between the amplitude of the equivalent white 
noise and the amplitude of the edge busyness for fixed threshold 
values. (It will be shown that their amplitudes are proportional 
to one another.) 

(iz) The relationship between the amplitude of the equivalent white 
noise and the threshold setting. (It will be demonstrated that 
the relationship is strongly dependent on the scene.) 

(zat) An expression for the net equivalent white noise when several 
unrelated distortions are introduced into the picture. (It will 
be shown that the power of the white noise that is equivalent to 
the combined distortions equals the sum of the powers of the 
white noises that, are equivalent to the individual components 
of the distortions.) 


II. EFFECT OF VARYING THE AMPLITUDE OF THE EDGE NOISE 


The amplitude of the white noise that was considered to be equivalent 
to edge busyness was measured for a number of amplitudes of the 
busyness. Typical results obtained by eight observers are plotted in 
Fig. 2. The ordinate of this graph is divided into four regions which 
correspond to subjective classifications of the noise: Objectionable, 
Annoying, Tolerable, and Undetectable. We are mostly interested in 
the region where the noise is tolerable; here the results are represented 
quite well by a straight line having unit slope. This indicates that the 
amplitude of the equivalent white noise is proportional to the amplitude 
of the edge noise. The constant of proportionality depends on the 
observer, but the variation is only +5 percent. 

The data in Fig. 2 were obtained for a threshold value of 0.07; that is, 
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Fig. 2—The white noise that is equivalent to various amplitudes of edge noise 
for one threshold. There were eight observers. 


the noise was added to the signal whenever it changed by more than 
7 percent of peak amplitude in a Nyquist Interval (0.5 us). Figure 3 
shows similar graphs for different threshold values. Notice that for 
zero threshold, the white noise equals the edge noise; this is a con- 
sequence of the method used for calibrating the system. As the threshold 
increases, the amplitude of the equivalent noise decreases but for each 
threshold, the amplitudes of the two noises, as indicated by the at- 
tenuator settings, are proportional to one another. 

We now define F(x) as the relative power of the white noise for a 
threshold z. It is the power of the white noise, V,,(z), that is equivalent 
to unit noise added to the picture at times when the video signal changes 
at a rate greater than 7, Ie., 


F@) =— 3° (1) 


IV. THE RELATIVE WHITE NOISE AS A FUNCTION OF THRESHOLD 


The relative white noise is plotted against threshold values in Fig. 4b 
for the scene shown in Fig. 4a. Also shown in Fig. 4b is the probability 
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Fig. 3—The white noise that is equivalent to edge noise for six threshold values. 
There was one observer. 
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Fig. 4a, b—Picture of a girl with a lamp: (z) The relative white noise plotted 
against threshold. (77) The probability that the threshold is exceeded. 
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that the slope of the signal exceeds the threshold value. The relative 
noise for this scene may be represented approximately by the expression 


F(x) = exp (—2/)) (2) 


where ) is a parameter. A similar exponential relationship is obtained 
for other scenes that contain a variety of detail. The parameter \ 
depends on the amount of detail in the scene: \ is small for flat scenes 
and larger for more highly detailed scenes. 

The relationship (2) does not apply to scenes made up of abrupt 
changes separated by regions of uniform brightness such as a page of 
graphics, as illustrated in Fig. 5. At very small threshold values F(z) 
falls rapidly with increasing threshold because it is determined largely 
by noise that is visible on the flat background, but when z exceeds 0.02, 
the noise 1s concentrated on the characters. Thereafter, further in- 
creases of threshold have less effect on F(z). Other interesting scenes 
with their graphs are shown in Vigs. 6 and 7. Notice that signals with 
slopes exceeding 0.25 of the maximum possible rarely occur. 


V. COMBINING NOISES FROM SEVERAL SOURCES 


In order to design and analyze systems, we need an expression for 
the visibility of noise that is added to the signal at times when its 
slope has a particular value. Our measurements describe the visibility 
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Fig. 5a, b—Picture of graphics: (¢) The relative white noise plotted against thresh- 
old. (it) The probability that the threshold is exceeded. 
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Fig. 6a, b—Picture of Karen: (7) The relative white noise plotted against threshold. 
(7) The probability that the threshold is exceeded. 
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Fig. 7a, b—Picture of stripes: (7) The relative white noise plotted against threshold. 
(17) The probability that the threshold is exceeded. 
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of noise added to the signal at times when the magnitude of its slope 
exceeds specified values: it is an accumulation of the required function. 
In order to obtain the required function we need to know how the 
visibility of a combination of noises depends on the visibility of its 
components. 

There is strong evidence that, for a combination of distortions, the 
net equivalent white noise has power equal to the sum of the powers 
of white noise that are separately equivalent to the component dis- 
tortions. This hypothesis has been checked by measuring the equivalent 
white noise for two uncorrelated distortions, separately and in com- 
bination. The results are shown in Fig. 8. The abscissa of this graph 
is the amplitude of the white noise that was equivalent to the least 
visible component; and the ordinate is the amplitude of the white noise 
that was equivalent to the combined distortions. Both values are 
expressed as a ratio of the white noise that was equivalent to the more 
visible component. The continuous curve is the result expected assuming 
addition of power. Evidently the hypothesis is a reasonable one, although 
there is a tendency for the more visible distortion to mask the less 
visible one to a greater degree than power addition predicts. 

The distortions used for these measurements were: 


(c) White noise added to the entire picture 
(it) White noise added to sections of the picture 
(212) Noise added when the slope of the signal exceeds a threshold 
value 
(tv) Noise added when the slope of the signal lay in a certain interval. 


NET NOISE IN dB 





-6 -5 —4 -3 -2 -1 0 
RATIO OF NOISE COMPONENTS IN dB 


Fig. 8—The white noise that is equivalent to a combination of two distortions 
plotted against the least visible distortion. Both noise amplitudes are expressed as a 
ratio of the more visible distortion. 


1504 THE BELL SYSTEM TECHNICAL JOURNAL, SEPTEMBER 1972 


When two distortions of the same kind were introduced into the picture, 
the result always corresponded with points on the continuous curve 
in Fig. 8. Such results are not plotted on the graph. 


VI. THE VISIBILITY OF NOISE AS A FUNCTION OF SIGNAL SLOPE 


The expression f(x) will be used to describe the visibility of unit 
noise power added to the signal at times when the magnitude of its 
slope is x. f(x) is defined such that, when rms noise a(x) is added to the 
signal, the net equivalent white noise has power given by 


N? = i " ONG) ae: (3) 


No attempt has been made to separate the visibility into two com- 
ponents dependent on the sign of x. If it is necessary to accommodate 
one scene, it will surely be necessary to also accommodate its mirror 
image; it therefore is unwise to make the noise dependent on the sign 
of x. The expression f(z) can be determined from the subjective measure- 
ments of F(x). F(x) is the power of the white noise that is equivalent 
to unit noise added to the signal at times when the magnitude of the 
slope exceeds x. Therefore, 





F(x) = | f(x) dev (4) 
or 
ja) = —2O = _P@). (5) 


Values of f(a) are plotted in Fig. 9a for the scenes shown in Figs. 
4 through 7. Notice that the visibility of noise tends to decrease as 
the slope of the signal increases. This decrease is caused by a com- 
bination of three effects: It is well known that the ability of the eye 
to resolve detail decreases near abrupt changes of brightness,**’” 
thus we expect the noise to be less visible for large values of x. Moreover, 
the fraction of the signal that changes rapidly is usually small, as 
Figs. 4 through 7 demonstrate; thus, for large values of zx, noise is 
added only to a small fraction of the scene. These effects are counter- 
acted to some extent by the fact that the eye tends to concentrate on 
regions of the scene that carry essential information’’, thus the bright- 
ness boundaries tend to be more significant than the overall probability 
of occurrence indicates. In general, there is not a simple relationship 
between the function, f(z) and the probability density of slopes, p(z). 
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Fig. 9a—The relative visibility of noise f(x) plotted against slope: (i) Girl with a 
lamp. (ii) Karen. (iii) Graphics. (iv) Stripes. 


These functions are compared in Fig. 9b for two scenes; a page full 
of printed characters and a page containing only one line of characters. 

Attempts to measure f(x) directly instead of deriving its value from 
F(x) have been unsatisfactory for practical reasons. It is much easier 
to detect the times when the slope of the signal exceeds a threshold 
value than it is to detect the times when it les in a narrow range of 
values. Viewers find that it is easier to measure F(x) than f(z). 


VII. OPTIMUM COMPANDING 


Knowing the visibility function f(x), it is possible to find distributions 
of the noise, a(x), that minimize the net equivalent noise given by 
eq. (3). Before discussing applications of eq. (3), we must emphasize 
that its derivation assumes that the noise in the picture can be expressed 
as a function, a(x), of the local slope of the input. This will be true only 
in special cases, because of storage properties of the integrator that 
are needed at the receiver. In general, integration will retain the average 
value of the noise, thus causing streaks to emanate from edges of the 
scene and persist into flat areas. The above discussion is valid only 
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Fig. 9>—A comparison of the function f(x) with the probability density of slopes 
p(x), for two scenes. 


when the average value of the noise added to cach edge is small. This 
is true for Differential Quantization and it can be true for systems 
of the type illustrated in Fig. 10a, provided the channel noise has only 
high-frequency components. 

The advantages of transmitting video signals as their first derivative 
are easily realized in digital systems in which the dominant source of 
noise is the quantization process. The average value of this noise can 
be kept small by placing a feedback loop around the quantizer, as 
in a conventional Differential Quantizer.’*"”” 

In order to set up equations, let us use the model of a Differential 
Quantizer that is shown in Fig. 10b. In this circuit, the impairments 
introduced by quantization are represented as an added noise,” n(Z). 
To represent companding of the quantizer levels, the differential signal, 
z, is compressed according to the function g(-) before noise is added, 
and expanded afterwards according to the function A(-), which is the 
inverse of g(-). 

The signal emerging from the expandor is given by 


y = h[g(z) + nl, (6) 


and when the distortion is small, n « g(x), it may be written as 


y = Alg(x)] + nh'[g(z)I, (7) 


DESIGNING DIFFERENTIAL QUANTIZERS 1507 
PRE—EMPHASIS COMPRESS NOISY EXPAND DE-EMPHASIS 


CHANNEL 
eee 


Fig. 10a—Companding the derivative of a signal. 


h’(-) and g’(-) are the first derivatives of h({-) and g(-). Because g(-) 
and A(-) are monotonic and inverses of one another 


Alg(@z)] = « and A’[g(a)]g’(x) = 1 (8) 
therefore, 
y(t) = x(t) + n(t)/g’(2). (9) 


Now, let u(t) and v(t) be the input and output signals of the system, 
then 


vi +7) = y® + of) (10) 
and 
z(t) = u(t) — v(t) (11) 
therefore 
vit + 7) = ul) + n@)/g'(2) (12) 
and. 
a(t) = u(t) — uli — 7) + n(t)/g'(a). (18) 
Of particular interest are conditions where 
x(t) & u(t) — ult — 7), (14) 


that is, where x approximates the slope of the input signal. Then the 
noise present in the output signal is a function of the slope of the input 
and the result in eq. (8) applies. For convenience, normalize the ex- 
pansions by assuming n(¢) has unit rms value, then the output noise is 


a(x) = 1/g’(2). (15) 





Fig. 10b—Model of a Differential Quantizer. 
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In practical systems, it is convenient to place constraints on the 
value of g(-). For example, to have fixed signal levels, let 


g(z) > 0, gO) =O and g(1) = 1. (16) 


We therefore require that 
a S70) fad : wy = a (17) 


With these properties, the noise given by eq. (3) can be minimized 
using variational calculus. (Similar functions are minimized in Refs. 
14 and 15.) The visibility of the net noise is a minimum when 


aol) = Fe) [ P@) dx (18) 
ie, 
a(x) & [f(a) 4. 


The optimum noise amplification factor a(x) is plotted in Vig. 11 for 
four scenes. Curve (i) shows that, for simple scenes, a compandor can 
be used to increase the noise continuously as the slope of the signal 
increases. For very busy scenes, curve (iv) shows that the edge noise 
can be made about four times larger than the flat-area noise. 


VIII. DIFFERENTIAL QUANTIZATION 


We have represented the quantization process as additive noise; 
this is an oversimplification, as is the approximation in eq. (14). How- 
ever, the conclusions do apply to practical Differential Quantizers. 
Their action is illustrated in Fig. 12, which shows a typical response to 
an input rising at a constant rate. The output voltage is seen to step 
up at each sample time by an amount equal to a quantization level. 
If the signal change is not matched exactly by an available level, then 
an appropriate sequence of levels is used to approximate it. In general, 
the noise introduced by this process is determined by the spacing of 
the levels that are used; and a well-designed quantizer will use the 
levels closest to the local slope of the signal. We will now assume that 
the quantization noise, a,(x), is proportional to the spacing of the 
quantization levels that are adjacent to the slope, x, of the signal. 
This function is plotted in Fig. 13 for an eight-level and a sixteen-level 
quantizer. These scales were selected by experimental trial and error; 
the eight-level quantizer'’’’* was optimized for pictures of faces and 
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NOISE AMPLIFICATION ap(x) 





0 0.05 0.10 0.15 0.20 0.25 
SLOPE OF THE SIGNAL, x 


Fig. 11—The optimum noise amplification factor a (x): (i) Girl with a lamp. (ii) 
Karen. (iii) Graphics. (iv) Stripes. 
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TIME ——» 


Fig. 12—The action of a Differential Quantizer. 
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0 0.05 0.10 0.15 0.20 0.25 0.30 
SLOPE OF THE SIGNAL, x 


Fig. 18—The quantization intervals for an 8-level and a 16-level quantizer. 


is inadequate for graphics. More emphasis was given to the reproduction 
of high-detail graphics in the design of the sixteen-level scale.’* The 
shape of the graphs in Fig. 13 corresponds quite well with the graphs 
in Fig. 11; the eight-level quantizer approximates curve (i) and the 
sixteen-level quantizer approximates curve (iv). 

The assumption that quantizing noise is proportional to the spacing 
of the quantization levels is true if a,(x) is taken to be the average noise 
for values of x in each quantization interval. In reality, the variation 
of quantization noise with slope is very complex and has much fine 
structure. Similarly, any attempt to match the curves in Fig. 11 by 
line segments, such as in Fig. 13, must necessarily ignore the finer 
variations of a(x). This is not a disadvantage because the finer variations 
are very dependent on picture content and gain settings. Indeed, it 
is advisable to match the quantizer to an average of a(x) in order that 
a range of scenes can be accommodated. If the intent is to make the 
spacing of the quantization levels a,(x) equal to the average value of 
a(x), there is an easy routine for designing quantizers. It is described 
in Appendix B. 


IX. VISIBILITY OF THE NET QUANTIZATION NOISE 


It is useful to have the ability to calculate an expression that describes 
the net impairment of the received picture; for this purpose we will 
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use the amplitude of an equivalent white noise. The equivalent noise 
for uniform quantization is easily calculated if overloading is neglected. 
The error at the decision time ranges with uniform probability between 
plus and minus half a step interval. Therefore, 


d 
N, = V3 ; (19) 
where d is the level spacing. The advantage of optimally companding 
the level of the Differential Quantizer is shown in Appendix A to be* 


| [ “pQ) al. (20) 


The predicted advantage over uniform quantization is 17 dB for the 
picture shown in Fig. 6. It is customary to separate the advantage into 
two parts. One part is a consequence of the fact that quantization 
levels need not be provided for slopes greater than certain maximum 
values. For example, no slopes exceed 0.25 in this scene. The second 
advantage is a consequence of tapering the levels; this is predicted to 
be 5 dB for the picture in Fig. 6, about 5 dB improvement has been 
realized by companding an eight-level quantizer.” A similar improve- 
ment has also been obtained by weighting the steps of a delta modulator.° 

A method for calculating the equivalent signal-to-noise ratio for 
particular quantizers is described in Appendix C. It predicts 50 dB 
signal-to-noise ratio for the sixteen-level quantizer and 42 dB for the 
eight-level quantizer when processing the picture in Fig. 6. Both results 
agree reasonably well with practical measurements. 


X. THE VISIBILITY OF NOISE IN DIFFERENT PICTURES 


White noise has been used as a standard for measuring impairments 
of a video signal and has proven to be very useful for evaluating dis- 
tortions of a given scene. Care is needed, however, in using it to compare 
the visibility of distortions of different scenes because the visibility 
of the white noise itself depends on the picture content.” Indeed’, it 
is very difficult to obtain numerical comparisons of the visibility of a 
given distortion in different pictures. Table I lists the main properties 
of the pictures used in this paper, together with an estimate of the 
visibility of white noise. For example, a given amplitude of white noise 
was judged to be 2.5 dB less annoying in the picture of Fig. 4 than it 
was in the graphics of Fig. 5. 


* This assumes that the signal-to-noise ratio lies within the range for which our 
results apply, ie, 40-53 dB. 
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TABLE I—MAaIn PROPERTIES OF PICTURES SHOWN IN 








Fias. 4-7. 
Scene Peak/RMS Visibility of Noise 
Girl with lamp (Fig. 4) 11 dB 0 
Graphies (Fig. 5) 11 dB —2.5 dB 
Karen (Fig. 6) 13 dB —1.5dB 
Stripes (Fig. 7) 12 dB —-1.5daB 





The visibility of noise is also dependent on the motion in the scene. 
We have estimated that white noise is 3 dB less visible on the face of 
a person talking normally than it is when the face is stationary. 

It has also been observed that, for the picture in Fig. 6, the edge 
busyness is about 2 dB less visible when the noise is synchronized to 
the scanning waveforms than when it is random. The advantage of 
synchronizing the noise is fully realized only in stationary scenes. 


XI. DISCUSSION AND CONCLUSIONS 


We have described a method for measuring and evaluating how the 
visibility of noise depends on the instantaneous rate of change of the 
signal. The visibility of quantization noise can be significantly reduced 
by making its amplitude depend on the slope of the signal. The optimum 
distribution is very dependent on picture content, therefore measure- 
ments have been presented for several different scenes. These scenes 
were chosen as being typical of important classes and their properties 
should be useful for characterizing system requirements. 

The measurements have been used to design quantization scales 
for Differential Coders and to predict their performances. The results 
agree very well with scales that were obtained experimentally and 
with measured signal-to-noise ratios. 

It would be useful to relate the visibility of noise, as measured for 
this work, to basic properties of scenes and to psychophysical properties 
of viewers. This has proven to be very difficult. For example, attempts 
to relate the visibility of noise to the probability density of slopes of 
the signal have been successful only for scenes of uniform texture, 
such as Fig. 7a. 

The work reported here relates to the rate of change of brightness 
in the horizontal direction. Comparable results could be expected 
for the vertical direction.’® It will be interesting to see similar measure- 
ments obtained for other types of pre-emphasis besides simple dif- 


DESIGNING DIFFERENTIAL QUANTIZERS 1513 


ferentiation. The experiments should also be repeated with other inputs 
such as color and audio signals, and frame-to-frame changes. 
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APPENDIX A 


Net Noise in an Optimally Companded System 


An optimally companded system is understood to be one in which 
the noise is distributed according to eq. (18), ie, 


ate) = 0) fF P@) de. (21) 
The net relative-noise power is obtained by substituting in eq. (3) to get 
No = [ a(x) f(a) dx = j [ f(z) ae (22) 

When there is no companding, a(x) = 1 the noise is given by 
N? = [ “acide (23) 


Therefore, the advantage of companding can be expressed as 


Adz | [ AQ) ae | (24) 


because the measurements have been calibrated so that 


[ GO d= ha 1. (25) 


Most viewers, when judging pictures, are conscious of two components 
in the noise: flat area noise which corresponds approximately to areas 
of the picture where x < 0.02, and ‘‘edge-noise’”’ which corresponds to 
z > 0.02. With optimum companding, the “edge-noise’”’ usually pre- 
dominates because 


[ aa i . f(x) de. (26) 
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APPENDIX B 


Determining Quantization Levels 


The objective is to design a quantizer for which the relative spacing 
of the levels a,(x) represents a given value of the function a (x) in 
eq. (10). An example is shown in Fig. 14. The function a (x) is plotted 
as a continuous curve and is extended to negative values of x by as- 
suming it is an even function. The rectangles superimposed on this 
curve have the following properties: First, their height is proportional 
to the average value of a(x) during the interval defined by the width 
of the rectangle. Second, the inclination of the diagonal is the same for 
all the rectangles; this makes their width proportional to their height. 

The quantization levels are set equal to the value of x that corre- 
sponds with the sides of the rectangles. Then the spacing of the levels 
is proportional, in turn, to the width of the rectangle, the height of the 
rectangle, and to the average value of a (x); thus satisfying our re- 
quirements. The number of levels in the quantizer is one more than the 
number of rectangles and this is determined by the inclination selected 
for the diagonals. Quantizers with given numbers of levels can be 
obtained by an iteration of the design method. The decision levels 
would be placed midway between the quantization levels.” 

A “worst case” design is required for systems that must accommodate 
a variety of scenes. For example, for small values of z, the levels would 
be chosen to accommodate scenes with properties given by curve (i) 
of Fig. 11, and for large values of z, the levels would be chosen to ac- 
commodate scenes with properties (iv). 








— 0.05 0 0.05 0.10 0.15 0.20 0.25 
SLOPE OF THE SIGNAL, xX 


Fig. 14—A method for selecting quantization levels. 
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APPENDIX C 


Determining the Visibility of Quantization Notse 


Having determined the values of the quantization levels it is useful 
to be able to calculate the visibility of net quantization noise. The net 
quantization noise power will be given by 


Ni = fo aayfee) de (27) 


where a,(x) is the rms quantization noise which is assumed to be pro- 
portional to the spacing of the levels. The constant of proportionality 
is given by eq. (11). If the nth quantization level is /n then 


_ (, = Uin-1)) 
Oq(X) ~~ 2V3 ) 


The constraints (16) can be interpreted as requiring a fixed number 
of quantization levels, which are numbered in ascending order according 
to their magnitude. Therefore 


| ee ea (28) 


Un 
Ni=s G — 1. f 1) ae | | (29) 
be ee 
Assuming that the companding is symmetrical, J, = —l, , then using 
eq. (4) and neglecting overloading we obtain 


1 


Dee 
Nae 12 


| arta ig) oe oS hep t) Fa) oF aon| (30) 


for an even number of levels 2N, and 


Ni= 35 bs (her — ILC») — FOO) | (31) 


for an odd number of levels (2N + 1). 

These expressions are easily evaluated when J, and F are known. 
Values of the quantization levels 1, are found by using the approximate 
method described in Appendix B because we do not have an easy 
method for selecting values of | that minimize N? . It has been found 
that levels, determined by the method described, give signal-to-noise 
ratios that are only 2 dB less than that predicted by eq. 14, that is, 
for optimum shaping of the noise using the continuous function ao(x) 
instead of the discontinuous a,(z). 

Important cases are ones for which the quantization levels are 
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sufficiently close that f(x) is effectively constant during each interval. 
Then the method shown in Fig. 14 for obtaining levels gives optimum 
companding a, = a) . Equation (30) shows that every quantization 
interval then contributes the same amount to the visible noise. 
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Modeling the Growth of Jumpers on 
the Main Distributing Frame 


By S. HALFIN, C. J. McCALLUM, JR., and M. SEGAL 


(Manuscript received April 18, 1972) 


The buildup of dead jumpers in the Alain Distributing Frame (DF) 
plays a central role in MDF problems; for example, a recent survey ranked 
dead jumpers as the number two problem on a list of the most frequently 
reported AIDF problems. In this paper, two models are proposed to 
quantify the buildup of both live and dead jumpers and to investigate the 
factors influencing the buildup. These models provide tools for the analysis 
and comparison of possible solutions to the buildup problem. 


I. INTRODUCTION 


The Main Distributing I'rame (IDF) in a central office building 
serves as the connecting point between the cable outside and the equip- 
ment inside. Conventional MDI*’s are iron or wooden structures with 
terminal strips mounted on each side. The two sides are termed vertical 
and hortzontal due to the manner in which the terminal strips are 
mounted. Cable pairs from subscribers’ stations are terminated on the 
vertical side of the MDF, while line and trunk equipments are wired 
to terminals on the horizontal side of the frame. In order to provide 
service to a subscriber, it is necessary to connect his cable pair to the 
proper line equipment. A frameman makes this cross-connection by 
manually stringing a wire, called a jumper, between the corresponding 
vertical and horizontal terminals. These jumpers, which can easily be 
100 feet long, are laid along horizontal shelves in the frame. As more 
and more Jumpers are added to the frame, these horizontal shelves tend 
to become crowded. 

When the service to a particular station is discontinued, it is necessary 
to manually disconnect and remove the corresponding jumper. Un- 
fortunately, disconnected jumpers are not always removed, which gives 
rise to dead Jumpers in addition to the lve jumpers. 

In the last few years, main frames which were originally conceived 
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over seventy years ago are becoming highly congested with live, and 
dead, jumpers.’ The problem of the buildup of jumpers on the shelves 
of the MDI has been discussed extensively and many solutions have 
been suggested and examined.”"® These proposals include preferential 
assignment of line equipments to cable pairs, multipling line equipments 
to several appearances on the horizontal side of the frame, spreading 
line equipments along the horizontal side, using several separated 
frames connected by tie cables, mechanization of the assignment 
records, and others. A natural way to evaluate such proposals is by 
means of a growth model for the buildup of jumpers in the MDF. 
Desirable features of such models are: 


(z) The models should enable the comparison of the amount of live 
and dead jumpers during a given period of time, for various 
conditions and methods of administration of the MDF. 

(17) The models should use as inputs either data which are easily 
available or which can be theoretically deduced from available 
data. Such data are: the physical description of the frame, the 
rate of connect and disconnect orders, distribution of lengths of 
the jumpers, which can be deduced from administrative practices 
(for instance, the method of assignment), ete. 


One would expect: such models to be able to clarify why some MDF’s 
suffer from dead jumper accumulation, while others, though apparently 
similar, do not. Also, the models should facilitate evaluation of the 
proposed solutions. 

In this paper, we investigate the jumper buildup in a frame over time, 
with emphasis on the dead jumper buildup. Two models are proposed 
to quantify this buildup. 

The first model is based on the assumption that the fraction of dead 
jumpers not removed from the frame has a certain functional form which 
depends on two presumably measurable parameters which are described 
later. Given values of these two parameters and the values of two other 
frame parameters which are easily determined, the model can be used 
to iteratively calculate the jumper buildup over time starting from some 
initial state. 

In the second model, it is assumed that a dead jumper remains on the 
shelf either because its assignment record is erroneous, or because the 
force needed to remove it is too large. Based on these assumptions, 
recursive equations are deduced for the quantity of jumpers on a shelf 
at any given time. 

In Section JI, the first model is developed and some numerical results 
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are given. Section III is concerned with the second model. Some con- 
cluding comments are given in Section IV. 


II. DEVELOPMENT OF MODEL 1 


2.1 The Model 


Let ¢ = 0, 1, --- denote discrete points in time. Time period ¢ will 
refer to the time period between time ¢ and time ¢ + 1. Let J(¢) be the 
total number of jumpers (both live and dead) in the frame at time ¢; 
equivalently, J(¢) is the number of jumpers in the frame at the beginning 
of time period ¢ before any service orders for that time period are 
processed. Let c(¢) denote the number of connects (i.e., service orders 
which require the installation of a jumper) that occur during time period 
t. Similarly, let d(t) denote the number of disconnects (i.e., service 
orders which require the disconnection and removal of a jumper) that 
occur during time period ¢. For ease in exposition, we assume that each 
service order calls for the connection or the disconnection of only one 
jumper. 

Using the terms introduced above, we can write J(é + 1) as 


Jé+1)=JSO +e) -bOdo, t= 0,1,--: (1) 


where @(t), defined as the fraction of the disconnects d(t) that are 
removed, accounts for the fact that disconnected jumpers are sometimes 
left in the frame. Clearly, 0 S @(t) S 1. 

It is convenient to introduce the substitution 6(¢) = 1 — a(t) into 
eq. (1) where a(t) represents the fraction of disconnects that are not 
removed. Note that 0 S a(t) S 1. This gives 


Jti+ 1) = JO + ce — dt) + alt) d® (2) 


where c(t) — d(t) represents the net gain in live jumpers during time 
period ¢. Similarly, the expression a(t) d(t) represents the increase in 
dead jumpers during time period ¢. It is clear that we favor the situation 
in which a(t) is near 0 (equivalently, 6(¢) is near 1). In such a case, 
the jumper buildup in time will correspond only to live jumpers added 
to the frame and the additional frame congestion brought on by the 
presence of dead jumpers will not be felt. 

We now simplify the model and remove the dependence of c(t) and 
d(t) on time. Let S be the number of service orders per time period and 
let G be the fraction of the service orders that result in net gain, 1.e., 


S = c(t) + d(t) (3) 
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» _ ef) — dh 


"et + dt) 


for allt = 0, 1, --- . We assume S > O. Note that in the usual case 
in which c(t) 2 d(t) we have 0 S G@ S 1. Using eq. (3) and eq. (4) to 


eliminate c(t) and d(¢) from eq. (2), we have 


Jit+ 1) = J() + GS + 40 — @Sa(?). (5) 


(4) 


2.2 The Parameter a(t) 


We now investigate the parameter a(t), the fraction disconnects not 
removed. Though S and G were taken to be time independent, this 
seems to be an unreasonable assumption to make for a(t). One might 
expect that a(t) will increase with time as the number of jumpers in 
the frame increases. That is, when there are few jumpers in the frame 
(i.e., J(¢) is small) we expect that most disconnects will actually be 
removed (i.e., a(t) is near 0) since the frame is certainly not congested 
and only the disinclination of a frameman to remove an occasional 
disconnect or an error in the records will lead to its being left in the 
frame. However, as J(¢) increases, the frame becomes congested due to 
jumper buildup, and we expect that, due to physical limitations, fewer 
and fewer of the disconnects will actually be removed (i.e., a(t) is near 1). 
The removal of a disconnect will become physically difficult if not 
impossible and could jeopardize the continued operation of nearby live 
jumpers. Thus, we define a(t) to be an increasing function of J(é) with 
a range of 0 to 1. Furthermore, we require that a(t) actually reaches 
the value 1, 1.e., at some time ¢* when J(t*) = K, we have a(t*) = 1. 
kK, the jumper removal congestion number, indicates the number of 
jumpers at which disconnects can no longer be removed. In practice, 
frames reaching this stage require complete replacement or extensive 
cleanup campaigns. Note that when disconnects are not being removed, 
the growth rate of J (é) is the full rate of the connects. The assumed form 


of a(t) is 
J) ¥ 
wy - [Bl sosx © 
{ Li, Ji)>K 


where y, a non-negative real number, indicates the rate at which a(t) 
approaches 1. To see the effect of the number y on a(t), refer to Fig. 1. 
For y = 1, a() increases linearly with J(t). For y = 2, a(t) starts more 
slowly with respect to J(é) but then accelerates. The extreme case in 
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J (t) 


Fig. 1—Graph of a(t) vs J(é). 


which y = © corresponds to the case in which a(t) = OforO0 StS #*; 
this is the situation in which all disconnects are removed until the live 
jumper buildup exceeds K jumpers. For y = 0, we have the extreme 
ease in which a(t) = 1 for all ¢, i1.e., no disconnects are ever removed. 

A number of factors influence the parameters K and y. Those factors 
determining K are primarily frame parameters determined by the 
hardware of the frame and its use. One such factor is the frame configu- 
ration which includes such items as the number of verticals on the 
frame, the number of shelves, the positioning of various terminal 
densities on the vertical and horizontal sides of the frame, and the 
presence of bends in the frame (if any). Another factor influencing K 
is the jumper assignment procedure which, in part, determines neatness 
of the frame. For example, one might intuitively expect that preferential 
assignment would smooth out the jumper congestion inherent in random 
assignment. Thus, a frame would be able to accommodate a great 
number of jumpers before dead jumpers could no longer be removed. 
Hence, K would be larger under preferential assignment. 

The factors determining y are generally related to the human en- 
vironment. One such factor is record quality. A frameman who attempts 
to remove a disconnect based on instructions derived from erroneous 
records may not be able to complete the task. In such a situation, the 
jumper in question will usually become a dead jumper. Other factors 
influencing y are the numbers of both framemen and frame foremen 
and their skills and motivation. 


1522 THE BELL SYSTEM TECHNICAL JOURNAL, SEPTEMBER 1972 


Note that the model assumes that all hardware features and adminis- 
trative procedures of a frame are tied together in the parameters K 
and y. This feature of the model allows us to treat radically different 
configurations, assignment procedures, etc. This abstraction is in con- 
trast to the work of Swanson* in which physical properties of the frame 
are considered explicitly. 

If we substitute the form of a(é) given in eq. (6) into eq. (5), we get 


S(t) |" 

WO) a GS oC ays| 9 ‘ Pat (7) 
\ J(t) + 30 + GS, LS at; 

Note that as ¢ passes t*, a(t) becomes 1, and the future growth in J (é) 
is independent of K and y. Given the parameters G, S, K, and y, eq. (7) 
allows us to iteratively calculate the jumper buildup over time starting 
from some initial value J(0). In numerical calculations, J(0) is usually 
taken to be 0 though it would also be reasonable to interpret J(0) 
as the number of jumpers in the frame at its cutover time, which might 
be a substantial positive quantity. Note that if J(0) is not equal to 0, 
then a(0) is a number greater than 0. If one desires to grow a(t) from 0, 
he may replace a(t) in eq. (6) (and then eq. (7)) by 


=O], ete (8) 


| 1, J(t) > K. 


2.3 Computer Program and Numerical Results 


t+) = 


a*(t) = 


A computer program was written to perform the iterative calculations 
of eq. (7). It is worth noting that for cases in which y is a non-negative 
integer, it is unnecessary to perform the iterative scheme of eq. (7) 
since J(¢) can be obtained directly by solving a nonlinear first order 
differential equation. 

The jumper buildup J(é) is now examined for typical values of the 
parameters S, G, y, and K. Let S = 300 service orders per day and let 
the gain fraction G = 0.05. Assume y = 2 which seems to be reasonable. 
Typical estimates for reasonable values of K tend to be in the range of 
fifty to one hundred percent greater than the terminal capacity of the 
frame. By terminal capacity, denoted C, we mean that design parameter 
of the frame which indicates the maximum number of live jumpers that 
the frame can accommodate. Tor example, C might equal the number of 
horizontal terminals of the frame. We shall take C = 80,000 jumpers | 
and K = 120,000 jumpers. Finally, let J(0) = 0. 
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The jumper buildup based on the above parameter values is graphed 
versus time (in days) in Fig. 2. The straight line in the graph indicates 
the buildup of live jumpers alone due to the gain GS. At the rate of 
GS = 15 jumpers per day, the terminal capacity C is exhausted after 
5334 days or approximately 20 years considering 270 working days 
per year. At this point of 5334 days, J(t) equals 446,033 jumpers with 
only 18 percent of these jumpers live. The point in the graph surrounded 
by a small box indicates the time t* = 3264. At this time, a(t) has 
become equal to 1, and J(¢) is growing linearly at the rate of the connects 
4(1 + G)S = 157.5 jumpers per day or 10.5 times the rate of the gain. 

Figure 3 considers various values of K for the fixed value of y = 2. 
In this graph, J(¢) is plotted on a logarithmic scale. Notice that J(¢) 
decreases, for fixed values of ¢, as we increase K. The line corresponding 
to K = o represents the case in which the jumper buildup is due to 
live jumpers alone through the gain GS. 

In Fig. 4, we consider the case of K = 120,000 with various values 
of y. Note that J(t) decreases, for fixed values of t, as y is increased. 
Note also that ¢* is an increasing function of y. We see that by controlling 
the factors that influence y so as to increase y, we postpone the time ¢* 
at which point we no longer remove dead jumpers. The case of y = © 
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Fig. 2—Jumpers vs time fora = 2 and K = 120,000. 
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corresponds to the ideal situation in which the jumper buildup is due 
only to live jumpers. 

It should be stressed that the functions J(f) given in Figs. 2, 3, and 4 
follow from our assumed choice of parameter values. To actually use 
the model and predict jumper growth in a frame, one must estimate 
parameter values for the particular frame in question. 


2.4 Use of the Model 


One suggested use of the model is to predict the “breathing time’’ 
that a frame has before total congestion sets in. Presumably there is 
some total congestion number (in jumpers) such that when J(¢) reaches 
this number, the frame can no longer operate but must shut down. 
This total congestion number should not be confused with K, the jumper 
removal congestion number. Given the values of the parameters and 
the current jumper buildup, the model can be used to predict the time 
at which J(t) equals the total congestion number. The time until this 
total congestion time is the “breathing time.” 

Any use of the model is predicated on our ability to make good 
estimates of the values of the model parameters S, G, K, and y. Data 
is readily available for the parameters S and G. However, K and y are 
not so easy to estimate. Though more work needs to be done on ways 
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Fig. 3—Jumpers vs time for y = 2. 
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Fig. 4—Jumpers vs time for K = 120,000. 


to properly estimate K and 7, it is sufficient to obtain an estimate of 
J(t) for two or more (nonzero) values of ¢t to determine these parameters. 

Numerous refinements of the model can be considered. One could 
take S and G as functions of time as they most certainly are. In addition, 
although the form of a(t) given in eq. (6) is a reasonable starting point, 
more realistic functions need to be considered. One such effort in this 
direction is the model developed in the next section. That model which 
is built on a slightly different approach introduces parameters that are 
easier to estimate. Both friction and the weight of jumpers in the frame 
are explicitly considered as well as the percentage of errors in the assign- 
ment records. 


III. DEVELOPMENT OF MODEL 2 


3.1 The Model 


We make the following assumptions: 


(c) The frameman removes a disconnected jumper if the following 
two conditions are satisfied: there is no error in the assignment 
record of the jumper, and the force needed to pull it out is not 
greater than Ry pounds. 
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(it) The force required to pull out a jumper is proportional to its 
length, and to the total weight of all the jumpers which were 
put in the shelf after it. 

(iz) The length of a jumper, the correctness of its record, and whether 
or not it will be disconnected at a given time are statistically 
independent events. 


Let us denote by X, the length of a live jumper which was installed 
at time ft, let F(x, t) = Pr(X, S x) be the distribution function of X, , 
and let X, denote the expected value of X,. Let B(é, r) be the expected 
total length of all jumpers at time ¢, which were installed in the interval 
[é — 7, t]. Let A(t, 7) be the expected number of live jumpers at time ¢ 
which were installed in the interval [tf — 7, ¢]. 

Next, let R(X, t, 7) be the force required to pull out a jumper of 
length X at time ¢, given that it is of age 7. According to assumption (22), 


R(X, t, 7) = R*(t, DX 


where R*(, 7) is the force required per unit of length. R*(t, 7) is arandom 
variable proportional to the total weight, and hence to the total length, 
of all the jumpers which are of age S +. 

To simplify matters, we make R*(é, 7) deterministic by replacing it 
by its expected value aB(t, 7) where a is a constant which depends on 
the weight per unit of length of a jumper, the dimensions of the shelf 
and the average friction of the jumpers (see Section 3.2). Note that the 
friction may be higher than that attributable to the theoretical friction 
coefficient of the jumper’s coating, because the jumpers consist of 
twisted pairs and do not generally lie in straight lines. Thus we assume: 








R(X, t, r) = aBlt, 7)X. (9) 
Let 
Ut, 0) = = me 5 (10) 


Thus X(t, 7) is the maximal length of a “pullable” jumper of age r. 

Let ¢ denote the probability of an incorrect record, and let H(t, 7) 
denote the average removed length of a disconnected jumper of age r 
at time t. H(t, r) is, according to assumption (777), the product of the 
expected proportion of such live jumpers having correct records and 
the expected length of such live jumpers, given that their length is 
bounded by X(t, r). Hence, 


Rt.r) 


Et, 7) = (-¢) | Pru coe eee (11) 
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We now introduce assumptions about the rates of connect and dis- 
connect orders. Let us denote the connect order rate at time t by X(t). 
A(t) may depend on A(t, t)-the total expected number of live jumpers 
on the shelf, which represents the size of the main frame, or on A(E, ty) 
(where say ¢, = 1 year)-the recent growth of the shelf, or on external 
parameters such as changes in the community which is served by the 
central office. 

Let u(t, 7) be the percent of disconnect orders for jumpers of age r 
at time t. This rate will be, in many cases, independent of t. If a jumper 
is of age r at time ¢, then the probability of a disconnect order for it in a 
small interval [¢, ¢ + At] is approximately p(t, 7) At. 

Now let At be a time increment small enough so that quadratic and 
higher powers of Aé can be ignored. We then obtain the following system 
of equations: 

[r/At] 


B(t + At, 7 + Ad) = Bil, r) — At >> E(t, iAdu(t, iAd 


i= 


-[A(t, At) — A(t, G — 1) Ad] + Ad) AEX, 
{r/At] 


A(t + At, r+ At) = A(t, 7) — At D7 alt, iAd) 


t= 


-[A(é, tad) — A(t, @ — 1) Ad] + AC) At (12) 


for all ¢, r such that t 2 + = O, with the initial conditions: 


B(t, 0) = ht ‘ 


A(t, 0) = 0 
BO, 0) — Aix a 
A(0, 0) = Ao 


where A, is the initial number of jumpers placed on the shelf. 
By letting At — 0, one obtains a system of partial differential-integral 
equations: 


But, 7) + Blt) = —] BU, Wut, AAG Ww + OX, 


AGs) PAS = [ Da LN: (13) 


Since £(¢, 7) is nonlinear in the unknown function B(é, r), it seems 
that an analytic treatment of these equations will not lead us far. On 
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the other hand, a numerical solution is readily obtained by using the 
recursive equation given in eq. (12). 


3.2 An EHxplicit Formula for a 


The constant a can be evaluated in the following way. Denote 

I—the length of the shelf. 

w—the width of the shelf. 

o—the weight of a pair of jumpers per unit of length. 

p—the number of pairs of jumpers in a unit of area in the cross- 
section of the shelf. 

6—the friction coefficient of the pairs of jumpers. (6 is the force 
needed to pull a pair per unit of vertical force which is applied 
on the pair.) 


We assume that the buildup of jumpers is homogeneous throughout 
the length and the width of the shelf, and that one jumper lies higher 
than another jumper if and only if it was installed later. Likewise we 
assume that there are wp pairs of jumpers in any horizontal layer 
at any cross section of the shelf. 

Let us consider a pair of jumpers of Jength X which was installed at 
time ¢ — +, and let the present time be ¢. Then the total vertical weight 
on the pair is 

x 4 
oB(t, 7) Ve 


Thus the force which is needed to pull it out is 


1 . 
wv p 





R(X, t, 7) = doBet, ys 


Comparing this formula with eq. (9) we obtain 
6a ; 
lwV p 


3.3 Distributions of Jumper Lengths 





a= 


One of the major factors determining the behavior of the main frame 
is the distribution of jumper lengths. One might expect that the shorter 
the jumpers, the slower is the growth of the mass of jumpers on a frame. 
There are two reasons for this phenomenon: (7) short jumpers have 
small mass, and (77) short jumpers are easier to pull out when discon- 
nected. The distribution of jumper lengths is determined, and hence 
can be controlled, by the method used for assigning vertical terminals 
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to horizontal terminals. Various administrative and planning methods 
have been proposed, such as preferential assignment, duplication of 
equipment terminals, etc. We believe that with the present model one 
can analyze the effect of any of these methods on the buildup of jumpers. 
To perform such an analysis one needs to determine the distributions 
F(a, t) for the various methods. We present here the distributions for 
the case of random assignment, and for the case of preferential assign- 
ment with leakage. 
Random assignment yields the following distribution: 


F(x, t) = 2 a(2l — 2) for OSaSl 


1 for l<z 


where / is the length of the shelf. Note that, in general, / will be a func- 
tion of t due to growth of the frame. 

By preferential assignment with leakage we mean that the MDF is 
divided into N zones with equal length //N. We assume that there are 
two types of connect orders: type 1, which can be connected within 
the home zone, and type 2, which, for various reasons, are randomly 
assigned over all zones. Let ¢ be the fraction of type 2 orders among all 
the orders. 0 S ¢ S 1. ¢ will be called the leakage coefficient. Our model 
is somewhat different from the model of Swanson.* There, the leakage 
is defined as the fraction of those orders which cannot be executed within 
one zone. Note that in our case, type 2 orders may still be executed in 
the home zone. If Swanson’s leakage coefficient is denoted by ¢, then 
we have 


N—1 
ees ae 





The distribution of jumper length for preferential assignment with 
N zones and leakage coefficient ¢ is: 





2 
ao (2p -«)+¢%@1—2) for 0<e55 
F(x, ) = 1-¢+¢5Q@l-—2) i re 
l N 
1 for Ll< xz. 


As before, the length / of the frame, the number of zones and the 
leakage coefficient may be functions of t, that is, 2 = l(é), N = N(é) and 
¢ = ¢(t). The case ¢ = 1 corresponds to the random assignment method. 
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Fig. 5—Growth curves for various values of the maximal pulling force Ro, and the 
fraction of erroneous records ¢, for random assignments. 


3.4 Computer Program and Numerical Results 


The model was programmed for the case of preferential assignment 
with leakage, where the length of the frame, the number of zones and 
leakage factor remain constant over time. As mentioned before, this 
includes random assignment as a special case. We made the following 
simplifying assumptions about the connect and disconnect processes: 


(zt) A(t) = AACE, 1), namely, the rate of connect orders is linear with 
respect to the total number of live jumpers. 

(21) p(t, 7) = wo , namely, the percent of disconnect orders is identical 
for all age groups of jumpers, and is independent of time. 


These assumptions may be changed by minor modifications of the 
program. 

The program solves the system of equations given in eq. (12) for a 
given time step size At. The user must specify \) At and py At. Note 
that At should not be too large, otherwise eq. (12) will not be valid. 
A good criterion for the size of At is that A(é, #) and Bit, t) should not 
be much different from A(t + At, t + At) and B(t + At, t + At) respec- 
tively (say by 1 percent at most). 
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Other input parameters are: Ry, a, l, ¢, &; Ao, N and finally n, the 
desired number of iterations, so that the equations are solved up to the 
time Tina. = n At. At the starting point for the iteration, = 0, the 
shelf contains 4, = A(0, 0) jumpers. 

The output gives, at various predetermined points in time, the number 
of live jumpers on the shelf, the total length of the live jumpers, the 
total length of all the jumpers, and the percentage of live jumpers in 
the mass on the shelf. 

Sample runs of the program were made with the following input 
parameters: 


At = 1 month 
No At = 0.063 
Mo At = 0.057 
n = 240 (20 years) 
A, = 1000 jumpers 
1 = 100 ft 
obit 
150,000 


One set of runs was made for the case § = 1 (random assignments) for 
values of ¢ ranging between 0 and 0.25, and values of R, ranging between 
0 to 67 lb. The results are illustrated in Figs. 5, 6, and 7. 
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Fig. 6—Percent of live jumpers after 10 years, as a function of the maximal pulling 
force Ry and the fraction of erroneous records ¢, for random assignments. 
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A second set of runs was made for fixed R, = 40 lb and fixed @ = 0.1 
with values of ¢ between 0 and 0.8 and N = 4, 6, 8. The results are 
illustrated in Vig. 8. The nominal data for this hypothetical shelf seems 
reasonably in line with some existing frames having approximately 
7.2 percent annual growth of live jumpers; / corresponds to an MDF 
with 150 verticals. 

The formula given previously in Section 3.2 can be used to evaluate 
the coefficient a. This formula requires knowledge of the friction coeffi- 
cient. Since references for the magnitude of this coefficient do not appear 
to be readily available, we resorted to the following alternative approach. 

By the same reasoning used to deduce eq. (9), we have 

ae 

aes: 
where # is the force needed to pull out the jumper, X is the length of 
the jumper, and B is the total length of all the jumpers which are higher 
on the shelf relative to this particular jumper. Thus, a can be estimated 
by pulling jumpers from the shelf, noting the required force, the length, 

and the depth of the jumper in the pile. 

Experiments of this type were conducted by Tengelsen:’ ‘Readings 
taken thus far range from 39 pounds for a 29 foot Jumper located deep 
in a pileup to a value of 9 pounds for a jumper approximately 100 feet 
long in the upper part of the pileup”. We estimate that ‘‘deep in a pileup” 
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Fig. 7—Percent of live jumpers after 20 years, as a function of the maximal 
pulling force Ry and the fraction of erroneous records ¢, for random assignments. 
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Fig. 8—Percent of live jumpers as a function of the leakage coefficient and the 
number of zones N, for preferential assignments with Ry = 40 lb and ¢ = 0.10. 


corresponds to B of 200,000 feet. Thus 


See ate ee Oy 
29 X 200,000 ~ 150,000 


Finally, in the same paper (Ref. 5) it is stated that a value of 40 pounds 
was recommended as a safe working limit. This gave us an approximate 
upper limit for Ry . 

From Fig. 5, one can deduce that for the range of parameters chosen, 
the accumulation of jumpers is more sensitive to ¢ than to R, . This fact 
is also reflected in Figs. 6 and 7. In I’ig. 6, we see that the percent of live 
jumpers on the shelf after 10 years is nearly independent of R, for 
R, > 30 lb. However, lig. 7 indicates that after 20 years, the magnitude 
of Ry is important. Figure 8 shows the influence of preferential assign- 
ment; we observe that dividing the frame into more than 4 zones seems 
to give no advantage when the leakage coefficient is fixed. However, 
one should be cautioned that, in practice, the leakage coefficient might 
be an increasing function of the number of zones. 


—_— 


Ib/ft?. 


IV. SUMMARY AND CONCLUSIONS 


Two models for the growth of jumper buildup on the MDI have been 
developed. Given the values of various model parameters obtainable 
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from data on frame characteristics and operations, these models permit 
the user to predict the jumper buildup on a frame. 

Many suggestions have been made for rehabilitating MDI ’s which 
are either jammed or deteriorating because of dead jumpers. Should the 
effort be in the direction of automating the records to reduce ¢, or in 
the direction of requiring less force to remove jumpers or enabling the 
framemen to pull harder? Does preferential assignment offer any advan- 
tages (for a particular MDF), and if yes, into how many zones should 
the MDF be divided? We believe that many of these questions can be 
answered with the aid of the models presented here. 
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Multidimensional Polynomial Algebra 
for Bubble Circuits 


By S. V. AHAMED 


(Manuscript received November 17, 1971) 


The paper expands the basic concepts of the coding theorists in the 
representation of data strings by algebraic polynonuals, and develops the 
representation of both time and location of individual binary positions in 
the same polynomial. Further, it advances a set of algebraic operations on 
such polynomials to correspond to the various subfunctions that are ac- 
complished in the actual domain circuits. The specific applications of the 
techniques proposed in this paper for the design and synthesis of such circuits 
ts presented in a companion paper. 


“Choose a set of symbols, endow them with certain properties and postulate 
certain relationships between them. Next, ... deduce further relationships between 
them .... We can apply this theory 7f we know the “exact physical significance’ of 
the symbols. ... The applied mathematician always has the problem of deciding 
what is the exact physical significance of the symbols. /f this is known, then at any 
stage in the theory we know the physical significance of our theorems. But the 
weakest link of physical significance is extremely fragile.’’ The original source of this 
principium is J. E. Kerrick in An Experimental Introduction to the Theory of Pro- 
bability, Belgisk Import Company, Copenhagen. It is also quoted in a slightly 
different form by F. M. Reza in An Jntroduction to Information Theory, McGraw- 
Hill Book Co., New York, 1961. 


I. INTRODUCTION 


Magnetic domains exist freely in thin platelets of orthoferrite crystals 
obtained by slicing them so that their crystalline axis is perpendicular 
to the surface of the platelets. Such domains are also present in very thin 
epitaxial garnet films (Fig. la) on suitable substrates. When the platelets 
or films are subjected to bias fields, these domains assume cylindrical 
shape and their diameter shrinks to microscopic sizes (Fig. 1b). Such 
domains (also called ‘‘bubbles’’) are stable under an appropriate bias 
field condition and they may be manipulated to perform’ storage, 
gating, looping and also certain elementary logic functions.” 

The domains are generally propagated from one location in the circuit 
to the next by subjecting them to the local bias field gradient. Basically 
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there are two methods of providing such a field gradient to propagate 
the bubbles. In the “‘field access propagation,’’* an alternating magnet- 
ization is imposed in a patterned soft magnetic overlay by an in-plane 
rotating magnetic field generated by a pair of coils carrying an alter- 
nating current. The coils completely surround the platelet with their 
axis in its plane. Two of the most commonly used overlay patterns are 
shown in Figs. 2a and b. During one cycle of the alternating current in 
the coils, the domains in the platelet move from one point in a pattern 
to the corresponding point in the adjoining pattern. This finite distance 
that the domain traverses during one cycle is defined as a ‘‘period.”’ In 





(a) 





(b) 


Fig. 1—(a) Magnetic domains as they are observed by Faraday effect in a typical 
epitaxial film 5 to 8 microns deep, deposited on Gadolinium-Gallium-Garnet (GGG) 
substrate 20 to 40 mils thick. Magnification 340. (b) Formation of “bubbles”? from 
magnetic domains at a bias field of 30 Oe in same material used in Fig. la. Magnifica- 
tion 340. 
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Fig. 2—(a) T-bar type of overlay used for field access propagation. (b) Y-bar type 
of overlay used for field access propagation. 


the “conductor propagation,” the local field gradient to move the 
bubbles is supplied by a current in a conductor (Tig. 3). A single phase 
current is generally used to periodically shift a bubble from one stable 
position to the next. Such stable positions are derived by a soft magnetic 
overlay also embedded on the platelet, and the periodicity of movement 
of the bubble from one position to the next depends on the frequency of 
the single phase excitation of the conductor. The distance which the 
bubble moves during one cycle of the single phase current is also defined 
as one “period.” 

Typical orthoferrites (YbFeO; , YF eO; , etc.) sustain 40 to 50 micron 
diameter bubbles, and the period is approximately 200 microns. Typical 
garnets (Er.Tb,Al,.,Fes.,0i2 and Gd,.,Tbo.7Fe,;0,.) can support 4 to 8 
micron diameter bubbles and the period is about 25 microns. The 
orthoferrites require about one micro-second to shift a bubble position 
by one period. The newer garnet materials also require about the same 
time, thus yielding a data rate of about one megacycle. It is customary to 
employ ‘‘bubble-no-bubble coding’’ with field access propagation and 
“lateral displacement coding” (LDC) with conductor propagation. In 
the former type of coding, the presence or absence of a bubble at an 
appropriate location denotes one or zero. In the latter type of coding, 
the bubble positions are coded as one or zero by laterally displacing 
them from one coding position to the other coding position (see Fig. 3). 

All the bits of information are propagated by one period in one clock 
cycle in the field access propagation. In conductor drive circuits with 
lateral displacement coding, the average velocity during propagation is 
generally limited to one finite value, even though information bits are 
sometimes held stationary. When one finite velocity of propagation is 
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Fig. 3—Conductor propagation of bubbles. 


assumed, it is possible to extend the capabilities of the conventional 
polynomial algebra used by coding theorists’’® to encompass the space- 
time relationships in bubble circuits. Bubble circuits perform in the time 
dimension and in the space dimension. The space dimension is further 
divided into subdimensions; the different sections (or elements) of the 
circuit which perform independently. Hence, any algebraic representa- 
tion should encompass the representation of time and the representation 
of space dimensions which constitute the circuit. It is proposed that the 
time dimension be associated with X and the space dimension be 
associated with Y in the algebra. 


II. DECOMPOSITION OF TIME, SPACE (CIRCUITS) AND FUNCTIONS 
Consider the subclassification of time and space dimensions as follows: 


The total time for a circuit to perform a function consists of a series of 
individual time intervals neccessary for the submodular functions. 
These individual time intervals can each be represented as a certain 
known number of clock cycles. Thus the unit of time is one clock cycle 
at the excitation frequency of the main field in field access propagation, 
or 1s one clock cycle at the drive circuit frequency in conductor 
propagation. 

The physical layout of the circuit can be classified into various sections 
(elements such as paths, loops, functional modules, ete.). Each element 
further consists of a certain predefined number of periods. This leads to 
the unit of physical (or spatial) dimension as one period corresponding 
to one pole pitch in the T-bar, Y-bar or chevron pattern in field access 
propagation, or to one pole pitch of the driving conductor in conductor 
propagation. 
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Next consider the space-time relation. In field access propagation, 
the speed of the bubbles is one period per clock cycle, and is influenced by 
the angular velocity of the main rotating field; and stationary bubbles 
are rarely encountered.* In conductor propagation, the frequency of the 
exciting circuit determines the velocity;* and stationary bubbles are 
commonly encountered. However, in every circuit, over a limited 
duration and within a preselected element of the circuit, one can express 
the space-time relation with absolute certainty. 

Finally, consider the overall algebraic representation. The entire 
circuit function is modeled by a series of algebraic operations, each one of 
which corresponds to a subfunction in the circuit. Each subfunction is 
carried out in the time dimension and in the space dimension. Hence, if 
we can resolve the function into its subfunctions, the circuit into its 
elements, and time into sets of clock cycles, and identify the individual 
subfunction with the circuit element and the appropriate set of clock 
cycles, then we can analyze and predict the functioning of the circuit 
with great accuracy. The representations of individual subfunctions by 
corresponding algebraic operations are developed in Section IV. 


III. REPRESENTATION OF TIME, LOCATIONS AND BINARY VALUES OF A BIT 
POSITION 


The origin of time may be chosen to be at any desired instant. How- 
ever, a certain amount of flexibility and ease of representation results if 
the origin of time is chosen to coincide with a definite function in the 
circuit. Generally, bubble circuits perform repetitive functions and it 
may be convenient to choose the origin of time at the start of a repetitive 
cycle. When circuits perform a wide variety of nonrepetitive functions, 
then the analysis should be attempted for each function independently 
to ascertain the correct operation of each one of the functions. In the 
algebraic analysis of bubble circuits, it is proposed that the exponent 
of X (the carrier of time dimension as introduced earlier) be used to 
represent the number of clock cycles that have elapsed between a 
prechosen origin of time and the instant under consideration. 

Further, it is proposed that the location of any given binary bubble 
position be represented by two components: (2) the element in which the 
binary position is presently located and (72) the exact period in that 


_ * The momentary variations of bubble speed at crossovers and compressors are 
ignored, and the entire distance is considered as one period. The effect of corners is 
dealt with separately. 
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element at which the binary position is currently present. This leads to 
designating a subscript to Y (the carrier of space dimension as intro- 
duced earlier) to indicate the element number and an exponent of Y to 
indicate the exact location within that element. 

Finally, it is proposed that the binary value of a bubble position at a 
given instant of time be denoted by a. The binary value changes as a 
bubble position passes through the known transition points in the circuit. 
However, at a given time, which is a certain number of clock cycles past a 
preselected origin of time (a known exponent of X), and at a given 
location (known subscript and known exponent of Y), the binary value of 
a bubble position is either known, or it can be determined with absolute 
certainty from other circuit considerations. 


3.1 Representation of an Isolated Bubble 


Examine a single binary bubble position (lig. 4a) the binary value of 
which is ‘a’ at an instant of time 7 clock cycles from a prechosen origin of 
time within the kth clement of the circuit located at the lth location. 
Then it may be represented as 


= ax'y! 

i 2 Base @ O41 
cc numa 
0 1 2 8 a 

u=axl vy, 
(a) 
LOCATION 
ALONG k th 
ELEMENT ———> 0, (oon See! 


eo OG 


VALUES ———~——®a,_; 3 An_p --- 


(b) 


Fig. 4—(a) Representation of a single bubble. See text for explanation of j,k, 
ie l. (b) Representation of a bubble stream j clock cycles after a prechosen origin of 
ime. 
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3.2 Representation of a Bubble Stream 


Consider a string of n binary bubble positions (lig. 4b). Let the binary 
values of these positions be do, a), -++ , @,-, during a clock eyele which is 7 
clock cycles from a prechosen origin or time. If these binary positions 
are located in the kth element at location J, , l,, --- , l,-1 , respectively, 
then the string of data can be represented as 


i=n-1 

w= QXi¥e +: a XY =X! DY a VE. (1) 
i1=0 

The sign of individual terms does not carry any significance. When the 

binary positions are adjacent to one another, then Jy, J,, +--+ , y-1 are 

consecutive numbers. 

Example 1: Consider four bubble positions the binary value of which 

18 M , @, , 22 anda; . If they occupy the 8, 2, 1, and O location of a sixth 

element after 28 cycles from a prechosen origin of time, then they may be 

represented as 


i=3 
u= X™*(aoYi + +++ a3Vg) = X” 2 ie Gaal (2) 


3.3 Lxplanation of the Algebraic Representation 


From the point of view of circuit analysis, the algebraic representation 
leads to the following propositions. 

(2) A series of bubble streams in a circuit are represented by a series of 
polynomials. 

(zz) Any one bubble stream is represented by a particular polynomial 
(the sum of individual terms). 

(iz) Each bit within a bubble stream is represented by a term (the 
product of components). 

(zv) The binary value of the bubble position is represented by a. 

(v) The number of clock cycles between a prechosen origin of time and 
the end of the cycle under consideration is the exponent of X. 

(vt) The location of the bit of information is represented in two 

sections: the element within the circuit (the subscript of Y) and the 
location within the element (the exponent of Y). 


3.4 Implications of Representation 
The implications of the prechosen representation are: 


(2) That a ‘‘snapshot” (i.e., a complete description of binary valucs of 
bubble positions and their respective locations) may be extracted from 
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the algebra after a predetermined interval of time from a prechosen 
origin of time (a known exponent of X). 

(iz) That there is a unique polynomial for every bubble stream in a 
circuit during any one cycle. 

(zz) That all the pertinent information regarding all binary bits of 
information in a stream is available within the algebraic polynomial 
representing it. 


IV. REPRESENTATION OF FUNCTIONS 


4.1 Generation of Bubble Streams 


Choose an origin of time at the end of the generation cycle of the first 
bubble position (i.e., as it is leaving the generator). If the generator is 
going to generate n binary positions the values of which are a , a, °°° , 
and a,-1 , then at the end of (n — 1) clock cycles, the binary string just 
generated may be written as 


t=n-1 

w= X™ May Y* + ++) aiY) = X™' YS a Y¥°? = (Ba) 
t=0 

at an instant when the a,,_, position is just leaving the generator. The 
alternate representation of the string after n cycles is 


t=n-1 


“w= xX” > aY°. (3b) 
1=0 


In (8a) and (8b), the exponents of Y indicate the locations along the 
bubble path which lead out of the generator. When the coefficients 
Mo, 4, °** , A, are consistently one, the action of an unconditional 
generator is represented (see Fig. 4 of Ref. 4, and Fig. 3 of Ref. 3 for 
T-bar and Y-bar configurations). 


4.2 Annihilation of Bubble Streams in Bubble-No-Bubble Coding 


The function of annihilation of any bubble stream u may be simply 
represented as 


u = 0. 


It is important to note that it is not identical to a polynomial in which 
all the binary bit positions are zero. When a conditional generator 
generates a string of binary bits which are all zero, it is still necessary to 
represent the binary string (3a) or (8b), since this string of data may 
interact with other strings at a later point in the circuit. 

The function of a conditional annihilator may be represented as the 
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action of a conditional gate (to be discussed later) and the action of an 
unconditional annihilator. 


4.2.1 Resetting of Streams in Lateral Displacement Coding 


In a functional sense, the function of resetting in Lateral Displacement 
Coding (LDC) is analogous to the function of annihilating in bubble-no- 
bubble coding. However, there is an important representational dif- 
ference. The location of a bubble position which is reset in LDC is still 
identifiable, whereas the bubble position in field access circuits loses its 
identity and location upon entering the annihilator. Upon resetting in 
LDC, the bubble positions enter a ‘‘dead interval’, or a sort of coma 
from the activity of the circuit, yet the algebra has to account for the 
elapsed time during which the bubble positions are in the reset state. 
Hence, a stream of n reset bubble positions in locations 0, 1, 2, --- , 
n — 1 may be represented in the usual way as 


i=n—-1 


i= xi »; a, (4) 


i=0 
where 7 represents the number of clock cycles from a prechosen origin of 
time, and where a , a, , --* , @,—-; invariably represent the reset status. 
Example 2: The representation of 3 reset bubble positions in location 9, 
5 and 1 after 20 clock cycles past a given time origin is 


1=2 
u = X*aY°? + a,Y + a,¥') = X” So a,Y°™. (4a) 
7=0 


4.3 Functions of One Stream Resulting in One Stream 


Let wu, be the initial polynomial about to undergo the function, and let 
u, be the resulting polynomial after m clock cycles. In a categorical sense, 
the operation may be represented as 


Uq = Balu); or, Fup) =a Ug ) 
where F may represent: temporary freezing, translation in the forward 


direction, translation in the reverse direction, looping, special looping 
used in memory operations for dynamic data allocation, ete. 


4.3.1 Temporary Freezing of Bubble Streams in their Locations 


This function is generally encountered in conductor pattern and rail 
propagation with lateral displacement coding. Consider an n bit data 
stream j clock cycles from a prechosen origin of time, and represented as 


t=n-1 


Up = ae be aY;, ? (5) 
i=0 
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where a and / are a; and 1; respectively. After freezing the movement of 
the stream for m clock cycles, we have the stream represented as 


t=n-1 


tg = X”-u, = X*™ DY aYi. (6) 
i=0 
If the bubble positions are all located at adjoining locations then 
t=n-1 
u, = Xi" YY aYi 4, (7) 
7=0 


where (/) is the location J) of the bubble position ap . 


4.3.2 Movement of Bubble Streams 


This is the most common bubble function and it is necessary to estab- 
lish a sense of directionality in the movement. If the bubble stream is 
moving so that a,-, bubble position of eq. (1) occupies position I, after a 
certain number of clock cycles, then a> would be the leading bubble. It is 
easier to work with a) as the leading bubble,* and this implies that 
lL >, >, +--+: > 1,-1. If the positions are in adjoining locations, then 
1;-, is 1; + 1. Consider a stream of bubbles represented by wu, in eq. (5), 
which has moved in the forward direction for m clock cycles, then the 
resulting polynomial wu, is 


t=n—-1 


Ug = X"Y"-u, = X*™ DY ai", (8) 
i=0 


where a and I as defined earlier are a; and 1; . If the data positions are in 
adjoining locations with dp» in (1), then 


t=n-1 
1 NO a (9) 
7=0 
Generally, the bubble stream crosses elemental boundaries when it 
moves. If this stream shifts from element k to element ¢ during the 
movement, and if their intersection is located at Yj; and Y? (Fig. 5), 
then 


t=n-1 


Ug = eS ay (10) 


t=0 


* This notation helps the circuit designer to comprehend the location and move- 
ment of the leading bubble position first, rather than comprehending the location and 
movement of the last bubble position first. 
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012 cee ag 0123 
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Z ONELEMENT k 
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Fig. 5—Forward movement with a shift from element k to element ¢. 


in general.* If the bubble positions are located in adjoining positions 
L=1,=)h—i= () —tand 


t=n—-1 
Uq = Xitm » Ge, (11) 
2=0 
Example 8: The representation of 7 bubble positions, the binary 
positions of which are a, through a, located at the 15th, 14th, --- , 9th 
locations of the 3rd element of a circuit at an instant of time, 36 clock 
cycles from a prechosen origin of time, is 


7=6 
Up = X*® Ya VP U. (10a) 
i=0 


Further, if this stream travels from element 3 to element 4, after travers- 
ing for 20 cycles with the boundary between elements 3 and 4 being 
located at Y3° and Y%, then the final polynomial is represented as 


i=6 

56 35-4 

Us, cms XxX iw a;Y3 ies 
t=0 


but 


35-7 35-7-20 15-7 
3 = V4 = Y, ) 


* An alternate way to visualize the crossing of boundaries is to write 


n—l 
tig = Xitm aY,'™, 
1-0 


and then replace Y,7+* by Y:2, yielding a = 1 + m — z, which leads to eq. (10). The 
prime indicates that the polynomial as such does not represent an observable stream 
but will do so after the next operation (s). 
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and thus 
7=6 
Foam (11a) 


10) 


i= 
The exponent of Y, becomes negative if the stream does not traverse at 
least 11 cycles when it is 7 bits long, or if the stream traverses for 20 
cycles but is 16 bits long. It is thus possible to have the stream segmented 
as u, and u, by an improper choice of the traversing time. Sometimes 
negative exponents may be carried for a few steps in the analysis without 
interpreting polynomials as streams during these steps. 
Examine a bubble stream represented by wu, in eq. (5) having moved 
in the backward direction for m clock cycles, then 
i=n-1 
igi OY a SX YS a Vi, 
i=0 
where a and / are a, , and 1; respectively, and dp is the last bubble in this 
ease. If the bubble stream moved from element & and entirely shifted 
in tth element, then 


i=n—- 


1 
w= Xi Ya, (12) 


where* the intersection of elements k and ¢ is located at Y? and Y?. If 
the bubble positions are in the adjoining locations, then 1; = lb — 7 
where Jy is the location of the bubble a) in u, before its transformation 
to u,. 


Example 4: The representation of a 4 bubble stream, the binary values of 
which are dp , @; , @2 and a; located in the 7, 5, 3, 1 periods of element 5 
at an instant 35 clock cycles after a prechosen origin of time, is 


t=3 
ti, ay ee (12a) 
i1=0 


After 9 cycles of backward movement, if the stream is in element 3 
with the boundary of elements 5 and 3 located at Y} and Y3° , then 


1=3 
i Dee (12b) 
t=0 


* An alternate way to visualize the crossing of the boundary is to write 
a t=n—1 
lig = XO DD a 
i-0 


and then replace Y;,* = Y.2+¢, when ais < 0, thus yielding z + a = z+ (1 — m) 
which leads to eq. (12). 
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4.3.3 Looping of Bubble Streams 


Looping is generally encountered in T-bar or Y-bar propagation when 
a fixed delay or storage is necessary for bubble streams without actually 
freezing their movement. 
Consider a bubble stream represented by u, in (5), which is circulating 
in a loop with g periods* for m clock cycles, then 
i=n-1 


p+m lim 
tS KS ye ys "se Gai 


7=0 


but in a loop (J + m = (1 + m) mod g) thus leading to 


t=n-1 


Uq = xe > 8 cea mod - (13) 
i=0 

where a and 1 are read as a; and /; respectively, and 1; = (lJ, — 2) if the 
bubble positions are in the adjoining location with a) as the leading 
bubble. 
Example 5: The representation of 26 bubble positions’ a), a, °** , Gas, 
occupying locations 28, 27, --- , 8 in a loop (element 5) at an instant 
1214 clock cycles from a prechosen origin of time is 


t=25 


ho Xe Sa, (13a) 
i=0 


After looping for 980 clock cycles, the final representation is 


t=25 t=25 


Ug = X2194 > a; Vo eocmenss Re bs Glue sy eye (13b) 


7=0 1=0 
4.3.4 Movement Around Corners 


Corners’ in T-bar and Y-bar circuits need special attention since 
bubble positions lose or gain a quarter period when a 90-degree turn is 
present. In most cases, their effect may be eliminated by considering a 
movement to span two or an even number of compensating 90-degree 
corners. However, when it is necessary to predict the movement of 


* The number of T-bar periods should be considered as the number of clock cycles 
to bring back the leading bubble to its original location in the loop. The orientation 
of the T-bars and the direction of rotation both play an important part in the deter- 
mination of the number. In any case, the value of g is the actual number of T-bars +1 
depending on the orientation of T-bars and the direction of rotation. 

tT One encounters this bubble stream in (39, 26) the shortened BCH encoder with 
ene domains constructed along the same principles as discussed in Section 3 of 

ef. 7. 

t When the effect of the corners in a loop is being considered, the value of k used 
should be (& mod 4), and the number of periods in the loop should be measured as 
indicated in See. 4.3.3. Else the value of & should be its real value. 
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streams which include k number of 90-degree unidirectional corners, 
then a polynomial wu, in eq. (5) becomes 


Ug = b 4a fae ‘Up (14) 


where the sign is determined by the direction of the turn with respect to 
the direction of rotation of the main driving field and the orientation of 
the T-bars. 

No special algebraic consideration is necessary in conductor pattern 
and rail propagation. 
Example 6: Consider the two bubble positions represented as 


Up = X” (ao Vi + a, Y$). (14a) 


If these positions traverse for 18 cycles, and encounter three unidirec- 
tional —90-degree turns in the path, then the final polynomial is 


iy DY ag ee Xa Va ai VY, (14b) 


The fractional exponent of Y indicates that the binary positions a) and 
a, are lagging 270-degrees behind their corresponding positions had 
there been no turns. If it is necessary to bring them back in phase, then 
an additional 2 clock cycle is required, and the polynomial u, would be 


= xe (GoYao + aia); (14c) 
4.3.5 Routing of Bubble Streams 


This function plays a critical role when it is necessary to transfer 
streams of binary information into a certain branch of a circuit at a node. 
The algebraic representation after this function is identical to the 
polynomial u, in eq. (10), when the polynomial crosses the boundary of 
one element k at Yj and enters another element ¢ at Y} . In the general 
case, when the bubble stream is channeled from an element p at Y%, 
and enters another element ¢ at Y{*” , u, may be represented as 


i=n-1 


Ua = xe ss i (15) 


t=0 


where a and / as defined earlier are a; and 1; respectively. 


4.3.6 Inverting the Binary Content of Bubble Streams 


Consider a polynomial 


t=n-1 


Uy = xe ys din¥; . (16) 


7=0 
If the bubble stream so represented has gone through an inverting gate 
during the following m cycles, then the resulting polynomial may be 
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written as 
. : 1 
Ug — Bis ee Ae Sasa (16a) 


where / denotes J; and 
Qig = ip. 


Example 7: Consider 3 data positions in an LDC circuit. If ao, a, a 
denote their status at an instant 49 clock cycles from a time origin, and 
are located in the seventh element at locations 17, 13 and 9 then, 


i=2 
liane baa s (16b) 
7=0 
If these binary positions pass through an inverting gate located at Y>° 
and Y;, then the bubble stream after 15 clock cycles is 


+=2 
ig = YOY. (16ce) 
7=0 


4.3.7 Opening and Closing Gaps in Bubble Streams 


This function is quite effectively used in dynamic data reallocation 
with T-bars. The bubble stream in the loop has two preferred paths, one 
for each direction of rotation. In one direction, the stream traverses 
(n + 1) periods in the loop, and in the other, it traverses (x) periods in 
the loop. (See Fig. 6). If after z clock cycles of clockwise movement 
of a bubble stream in the loop the direction is reversed for one clock 
cycle, then the data bit at the nth (lig. 6) location is at the Oth location, 
and the data position which was at ((¢g — 1) mod n) location is at the nth 
location. Now z clock cycles of anticlockwise rotation would have 
effectively included the data position at the nth bit at the (2 mod (n + 1)) 
location, and all the remaining positions one location behind their 
original locations. A converse process takes place for excluding a data 
bit position in the stream. 


The algebraic equivalent of this function may be represented as 
follows: 


Let the numbering of the locations in the loop be in the direction 
(anticlockwise in Fig. 6) which permits the stream to traverse n periods 
in the loop. Let the contents of the loop (without the nth position) be 
represented as 


t=n—-1 
(i= I a (17) 
i=0 
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Fig. 6—Schematic representation of the dynamic data reallocation function. 


The polynomial u* undergoes a shift in clockwise direction for m (with 
m <n — 1) clock cycles. The resulting polynomial uw, becomes 


t=n-1 


Ud = D Glis Ps any rood (18) 


t=0 


As a distinct step from the above representation, let 
Up KX Do a V" (19) 
1=0 


represent the contents of the loop including the nth location before a 
series of counter clockwise shifts for m clock cycles. The resulting poly- 
nomial is 


n 
Uy = xm Ne ee Some mod Ae Be (20) 
0 


Example 8: Consider a string of bubbles the binary values of which 
are A, 4, °*: , @ at an instant 92 clock cycles from a time origin 
located at 0, 1,2, --- , 7in a memory loop (I'ig. 6) with 255 (i.e., nm — 1) 


* There is no need for a subscript for Y since only one element (the loop) is being 
considered. 
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periods in the clockwise direction and 256 periods in the anticlockwise 
direction. This can be represented as 


i=7 
Up = X” DYoaY'. (17a) 

i=0 
Let the 256th (i.e., nth) position contain the bubble position a, at the 


same instant of time. Four cycles of rotation of the field causing a 
clockwise circulation of binary positions yields 


t=3 1=7 
Ug = aps a; yet ae > a, yd. (18a) 
i=0 t=4 


In the first four terms in the polynomial (38), (¢ — 4) mod 256 would 
correspond to 252, 253, 254 and 255 for 7 ranging from 0 to 3. 

Now if the field is rotated for 5 clock cycles in the opposite direction, 
resulting in an anticlockwise shift of bubble positions, then the resulting 
polynomial is 


i=3 i=7 
= or Py a, mod 257 Yo ie ; 
U > a +a + dX a 
or 
1=3 i=7 
“= ra> a:Y'+a,Y'+ >> ave. (20a) 
1=0 i=4 


In effect it is seen that the bubble position after the nth location has been 
inserted at the 4th location in the bubble stream. 
4.4 Functions of One Stream Resulting in Two or More Streams 


Let u, be the initial polynomial about to undergo the function resulting 
in two polynomials uw, and u, . In general, functions of this type may be 
represented as 


Ug oF Uu, = Ha Us) or, Fits) ae Ug - U, ) 


where F* represents duplication or addressing sections of the initial 

polynomial u, into one or the other branch of a circuit. A finite number 

of clock cycles (m) are allowed during which the operation takes place. 
4.4.1 Duplication and Replication of Bubble Streams 


Duplicators may require a certain finite number of clock cycles to 


* These functions are used extensively in the companion paper® dealing with the 
applications of the algebra. 
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Fig. 7—The bubble stream uw, entering the duplicator D. 


operate. In the T-bar propagation, one clock cycle is necessary to dupli- 
cate. The effect may easily be included in the algebra by considering the 
duplicator as an extra period of element p with u, (see eq. 5), and the 
zero location of elements q and r in which u, and u, will be positioned 
after m clock cycles. (See Fig. 7). 

Consider a string of binary data represented as u, in eq. (5) to be 
duplicated during the m clock cycles yielding u, and u, represented as 


t=n-1 


ip. Se ea RO Say (21) 
7=0 
i=n-1 

te = XY, = XM ayers, (22) 
i=0 


where a and / are a; and I; respectively, the duplicator is located after z 
periods in k, and d is the number of clock cycles to accomplish the 
duplication. It is to be noted that when the bubbles occupy adjoining 
locations, d cannot exceed 1 for satisfactory duplication. Under such 
conditions* 


t=n-1 


Uq = yite >> a et (23) 
t=0 
and a similar expression for u, , where (/) denotes Jy . 


* An alternate way to visualize this transformation is to consider that the dupli- 
cator is located at the intersection of Y,?t!, Y,° and Y,° which leads to Y,* = 
Yr tv = Y,e-G@t), when a > z + 1, thus yielding 


t=n—-l 
Ug = Xitm™ FS aYzit)-itm, 
t-0 


where (l) is the location lp of the first bubble position a» prior to duplication. This 
leads to u, in eq. (23). 
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Replication leads to a series of resulting polynomials u, , u,, --- , ete., 
but the algebraic treatment is exactly the same. 


4.4.2 Addressing Sections of a Stream into Different Branches of a Circuit 


Consider a polynomial u, in element p which will approach a gate. 
The gate operating for g clock cycles will address the first n’ bubble 
positions in qg and the rest into 7. This function can be visualized as the 
effect of two independent translatory functions: (see Section 4.3.1) (2) 
the first n’ positions are translated in the forward direction and change 
elements from p to q, and (iz) the remaining positions are translated in 
the forward direction and change element from p to 7. Algebraically this 
can be expressed as follows: 


t=n—1 t=n'-1 t=n—-1 


w= Xi Voayiax' S ayi'+ xX’ ¥ av? (24) 


124 =3 can! 
with 

Ug = xi*™ pa aver (25) 
and 

u, = Xi" ayes (26) 


t=n’ 


where a and / are a; and 1; respectively (see Section 5.3.1). The gate 
should be located at Y3 , Y? and Y? . It can be seen that the first bubble 
position does not reach the gate till (eg — (Jo + 1)) clock cycle, and to 
divert the first n’ positions the gate should be operating to divert into 
element q for exactly (l, — 1, + 1) clock eycle, leading to the design 
detail that the gate should divert into q for g clock cycles where g = 
(2 — l, + 1). If the bubble positions are in the adjoining location, then 
the gate has to operate diverting into q for n’ clock cycles starting after 
(2 — (ly + 1)) clock cycles. Further, it has to act for the next (n — n’) 
clock cycles to divert the bubble position into r. 

When a gate addresses various sections of a data stream into more 
than two elements, the algebraic representation is similar. Such a con- 
dition exists if m is not chosen large enough in the previous case, and 
then there will be two bubble strings u, and u, together with a section 
of u, , which has not been processed by the gate. 

Example 9: Consider a data stream* in element 0 of a circuit. Eight data 


_ * Such a data stream is encountered in general rate change circuits represented 
in Fig. 1 and la of Ref. 9. 
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bits a) through a; , 12 clock cycles after a prechosen origin of time, 
occupy location Y;? through Y> . The polynomial representing the stream 
is 


i=7 
GS Saya (24a) 
7=0 


If the data stream passes through a gate at Y*? which diverts the first 4 
bits into element Y, and the last four into Y, , then after 14 clock cycles 
the resulting polynomial u, in 1 and wu, in 2 may be represented as 


i=3 
Ug = X* Da Vi, (25a) 
7=0 
and 
i=3 
KO ave ® (26a) 
i=0 


4.5 Functions of Two or More Streams Resulting in One Stream 


Consider two streams u, and u, (Fig. 8) interacting to yield one 
stream u, at a gate @), where @) may denote the function of logical 
gating, combining, etc. This function is different from the previous 
functions, since the individual binary values a; in the polynomials are 
likely to be changed by this function. Let 


t=n—-1 


uy DT as, (27) 
= 
and 
t=n’-1 
Ug am Gs > G2V (28) 
i=0 
ELEMENT p 
Un—1,p) {(0,p) 
U 
X+— p -4 
an-1 4op | 


OL—--« 


| Z’ ELEMENT r 


x — (x) (0) --- == : 


1 Cet) 0.4) 
ELEMENT q 


Fig. 8—Gating or merging of u, and u, to yield u,. 
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where [{zp) and liq) represent l;, , and 1;, respectively. The number of 
locations including the first and the last binary position of u, is (lop — 
lin-1.p) + 1); (i.e., n if adjoining locations). It can be seen (ig. 8) that 
the length of the resulting bubble stream, wu, , is the largest number n’’ 
of the following four numbers* 


(2) lop — Un-1.my + 1, 
(it) Vg dha ay Seal, 
(ttt) (lop — Utn-1.») 1+ 2 — lop — @ = boa) 
= (@— la-am) — (@ — log) + 1, 
or 
(iv) (log = en -1,09) 1 + 2’ — log — @% — boy) 
= (2! = devi.) — (@ — lop) + 1, 


where z and 2’ denote the number or periods along elements p, and q at 
which the gate © is located (see Tig. 8). The polynomial u, can be 
written as 


ioe > ag, (29) 


where /{i7) represents 1,;, respectively, and 


Air = Gin (Q) Gia (30) 
lL,~= lo fm — Zz; (81a) 

or 
lL, =l,, +m — 2’. (31b) 


The subscripts for a and / should be chosen with adequate care to 
consider only the interacting binary positions in streams u, and u, that 
pass through the gate simultaneously. It is important to note that 
each term in wu, results from a term in u, , or in u, , or from terms in both. 
When there is no term in one polynomial (u, or u,) corresponding to a 
particular term in the other (u, or u,), then the appropriate equation 
(81a or b) for the exponent of Y should be chosen. When there is a term 
in one polynomial (u, or u,), corresponding to a given term in (u, or u,), 


* (2) or (22) indicates the length wu, or wz with the longer stream uz, or u, completely 
overlapping w, or wp respectively. (277) or (7v) indicate partial overlap, with wg or up 
being nearer the gate. 
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then the location calculated from 1,, or 1;, will yield the same result for 
bets 

Example 10: Consider a bubble stream u, with four binary bits ao, , 
Gy, A, and a;,, after 41 clock cycles from an origin of time located at 
3, 2, 1 and 0 of element 5 and represented as 


7=3 
=X” Stns (27a) 
i= 


Consider a next stream u, with seven bubble positions uw, through us, ; 
also 41 clock cycles from the prechosen origin of time occupying 10, 9, 


8, --- , 4th locations of element 6, and represented as 
i=6 
je = KF a YS. (28a) 
7=0 


If u, and uw, approach an exclusive-or gate @ located at Y?, Y 4’, and 
Y?, then the polynomial uw, after 12 cycles can be derived as follows: 
Qoq » Gig pass through the gate during the first three cycles while u, is 
still approaching the gate. The binary positions a., through a;, interact 
with ao, through a3, in the gate @® for the next four cycles. a,,. passes 
through the gate and the gating is now complete. During the last five 
cycles, the bubble stream in element Y; moves away from the gate thus 
leading to the resulting polynomial 


i=6 
tae” Do Gave. (29a) 
i=0 


where do, , Qi, and ds, aYe Ao, , Ai, and a,, respectively, and 
Qi, = Ay Ba;, ¢ = 2 through 5). (30a) 


This example corresponds to the case (2) in Section 4.5. 


4.6 Functions of Two or More Streams Resulting in Two or More Streams 


This function, though rarely encountered in normal bubble circuits, 
can still be conveniently represented as an integral procedure of many 
subfunctions in which two or more streams result in one stream. If g, , 
Je, 93 *** are individual functions yielding streams 1, 2, 3 --- , etc., 
then the algebraic representation of Section 4.5 may be extended to 
represent streams 1, 2, 3 --: , etc. One such example is presented in 
Section 2 of Ref. 8. 


V. OVERALL BUBBLE CIRCUIT FUNCTIONS 


We have a set of mathematical tools to predict the binary values and 
locations of individual bubble positions as the binary streams undergo 
different submodular functions within the circuit. The interval of time 
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chosen for these submodular functions in Section IV is ‘m’ clock cycles. 
To study the overall circuit function, the function is divided into a series 
of submodular functions (7, , F, , F'3 -- +); the circuit which performs the 
functions is divided into a series of elements 


(12,3 oer g hy ke dy see), 


and the time necessary to accomplish the function is divided into a series 
of (m, , Mz ,m3 --+ , ete.) clock cycles. 


5.1 Subdivision of Circuit nto Elements 


After isolating the subfunctions within the overall circuit function, 
the elements that accomplish these subfunctions may be identified. A 
series of fine functional subdivisions may be necessary to identify the 
particular circuit elements. 

Some of the examples of the elements are storage paths, transmission 
paths, loops, ete. Sometimes a particular element (‘) of a circuit designed 
is very short, and it cannot accommodate the entire string of data. When 
it is still desired to study the contents of the n bit binary string of the 
polynomial u, in that element, then it is possible to fictitiously extend the 
element to just accommodate the n data bits. With the observer located 
at a preselected period (say Y;) in the element, the binary values of data 
which flow past this location would still be the values of ao , @ , dz +++ Gy 
in the calculated polynomial wu, . The instant of incidence of the leading 
bubble a) would be (n — b) clock cycles prior to its value as predicted by 
the exponent of X(i.e., jo) in the polynomial u, . 


5.2 Subdivision of Time 


Subfunctions are accomplished by elements within specified intervals 
of time. The interval of time for a specific subfunction is almost entirely 
determined by the circuit parameters and the clock frequency. Any one 
particular interval of time may, however, be conveniently expressed as a 
certain number of clock cycles. 

These different values of clock cycles m,, m2, m3, --: , ete., are 
necessary to calculate the polynomials wu, , U2, v3, °-: , etc., and they in 
turn uniquely define the values, intervals and positions of binary bits in 
the circuit. Generally, the subfunctions may proceed in series or in 
parallel, and different bubble streams may simultaneously undergo 
different functions in different elements of the circuit. However, the 
complete function of a circuit starts and finishes at an instant of time. 
Hence, whenever subfunctions are accomplished in parallel, it should be 
realized that the summation of m,; does not equal the total number of 
clock cycles necessary for the complete circuit function. Each sub- 
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function proceeding simultaneously should be modelled independently 
by its corresponding algebraic operation. 


VI. CONCLUSIONS 


The operation of bubble circuits may be effectively analyzed by multi- 
dimensional polynomial algebra without actually constructing the 
circuits. The location of all the data positions in the circuit can be 
accurately predicted at any preselected instant of time during the 
operation of a circuit by this technique. When all the circuit parameters 
are not accurately known, the analysis helps in the calculation of some 
of the circuit parameters. Further, it helps to algebraically check the 
validity and effectiveness of a conceived circuit in the performance of 
specified functions. 

The algebra may also be used for circuits that do not perform instantly 
but need a certain predetermined duration for movement, duplication, 
gating, etc. Several technologies (magnetic domain, charge coupled and 
charge transfer technologies) presently being developed in the Bell 
System fall into this category. 
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Applications of Multidimensional 


Polynomial Algebra to Bubble Circuits 


By 8. V. AHAMED 
(Manuscript received November 17, 1971) 


The principles of Multidimensional Polynomial Algebra developed 
in @ companion paper are applied to two T-bar circuits with bubble-no- 
bubble coding and one double rail circuit with lateral displacement coding. 
The object of this paper is to indicate the flexibility of the algebra in its 
use with real circuits and to emphasize the potential of the algebra as a 
design tool for bubble circutts. 


I. INTRODUCTION 


The operation of bubble circuits depends on the accurate functioning 
of individual elements such as channeling gates, logic gates, generators, 
etc., at the critical instants of time. When the circuit becomes com- 
plicated, it is not easy to comprehend a multiplicity of functions and 
determine accurately the instant and duration of operations of these 
critical elements. Further, a bubble circuit cannot be easily altered 
like a prototype experimental electrical circuit; and it becomes necessary 
to ascertain the proper functioning of the designed bubble circuit 
prior to its actual construction. 

In this paper the technique developed in Ref. 1 is applied to three 
circuits, and the operation of each circuit is predicted. In the second 
example, the design parameters such as the operating time of the gates, 
their duration, and the overall timing for the circuit, are derived step 
by step as the algebra progresses. 


II. APPLICATION OF ALGEBRA TO A T-BAR STATION SCAN MEMORY 


2.1 The Principle of Operation 


Consider a hundred stations with 25 lines each. The status of each 
line is stored in a loop* with 2500 periods as shown in Fig. 1. A controlled 


* This configuration of the station scan memory was supplied by A. J. Perneski 
and R. M. Smith. 
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Fig. 1—Block diagram of a station scan memory. 


generator G codes the status of each line (active or inactive) by gen- 
erating a bubble or no bubble during the clock cycles through which 
the line is being scanned. Thus each line is scanned every 2500 clock 
cycles. The operation of the gate g is based on the repulsion between 
bubbles (if any) arriving through elements 1 and 2. One bubble in 1 
or 2 passes through the gate to elements 3 or 4 respectively. Two bubbles 
in 1 and 2 are diverted into the annihilator A and into the element 5. 
The gate g thus stores the most recent status of each line in the loop 
through either element 5 or element 3. Two sensors A and I in paths 
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3 and 4 indicate whether the line has just become active or just become 
inactive. 

The six elements in the circuit are each marked serially so that the 
bubble streams can be observed in these elements. The gate g does 
not have an independent status as an element since the streams can 
only pass through the gate and they cannot be observed within it. 


2.2 The Algebra of the Circuit 


The gate g permits the interaction of the two bubble streams in 
elements 1 and 2, yielding streams in elements 3, 4 and 5. The truth 
table for the operation is: 





Inputs Outputs 
12 3845 
00 000 
01 010 
10 100 
11 001 


Let the generator G be located 20 periods behind the gate g (including 
its own location) on element 1. A binary position just generated by G 
would interact with location (/,), 20 periods behind the gate on the 
element 2. Further, let the elements 3 and 5 contain 30 periods each 
and let the sensors I and A be located 8 periods from g on elements 
4 and 3, respectively. These circuit conditions may be algebraically 
represented as: 


g located at Y?°, Y°*?°, Y3,, Yz and Y?} 
P located at Y3° and Y2° 
I and A located at Y{ and Y; , respectively. 


The circuit operates on a repetitive basis, and it is possible to choose 
an origin of time at the end of the generation cycle of any one bubble 
position at G. The origin also corresponds to the end of a coding cycle 
for a particular line L. The algebra may be carried out with any number 
of bubble positions in the polynomials. If four positions are chosen 
to illustrate the use of the algebra, then the four bit string after three 
cycles from this prechosen origin of time may be written as: 


Uy = X*(do1 Yy a are ae oF ce a = a3, Y4). (1) 


At this instant of consideration, the corresponding string of four bubble 
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positions in element 2 may be written as: 
U2 = X* (aoe Yo + ay Year + 22 Lee + Ase Yo") (2) 


It will take 20 clock cycles for the generated binary position a3; 
to be completely processed by the gate g. Thus with m, = 20, uf and ui 
(see footnote in Sec. 4.3.2 of Ref. 1) can be written as: 

u= X”*(do1 YO + ce ia + der ei + tsi Y i), (3) 


and 
23 lo +20 lot19 Lo t18 lo +17 
uw = X (doo Y; > + Qis Y: = + d2Y2 . + 32 Y; ‘“ y (4) 


The last positions a3, and a3, would be processed by the gate at the end 
of the 23rd cycle. Now the conversion of ui and uf to us , us and us is 
feasible by the truth tables for inputs and outputs, and by the spatial 
conversions: 


VP VE VE or Ye: 
20 
ga ka. SOR es 


Hence, if ao; , Qi, , Go, and as, are 0, 1, O and 1; and doe , da, , Ae. and 
Gz; are 1, 1, 1 and 0, then: 


Us = X”* (dos Y3 + ay; ae + Aes Y3 + 33 a) (5a) 
Us = XP (an Me af ayy, + Ao Y; + As4 YD (5b) 
Us = XG; + ae oe + re ae “5 Gans) (5¢) 


where G3 , Qi3 , Ge3 and az33 are 1, 0, 1 and 0; doa , Ais , Gog and G34 are 
0, 0, O and 1; and aps , dis , G25 and az; are 0, 1, 0 and 0, respectively. 

When u; and u, are multiplied by X°Y’, the exponents of Y; and Y; 
are both 8. The sensors A and I can read the status during the 28th 
(i.e., 23 + 5) clock cycle from the prechosen origin of time. In this 
case, the origin of time corresponds to the start of the coding cycle 
for the particular line the status of which has just been read, and the 
28 cycles indicate the delay between coding the status at G, and reading 
the change in status at A or I of any particular line L. 

Now consider the merger of the streams 3 and 5 at P, ie., at Y3° 
and Y;°. The merger is complete when the exponent of Y reaches 30 
which leads to 


U >= xX”?°Y*. (us -+ Us); 
or 


Us = X”*Aog Ye + ai6 ¥; + oe Ys + Avg Ve (6) 
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With og , Gig , 2g and dz, being 0, 1, O and 1. The position ao, reaches 
the merger point Y§ during the 50th cycle, a,, during the 51st. cycle 
and so on. 
If the Y% is 2447 (i.e., 2500-53) periods* from Y3°, then 
Uo = eae y744? ay 
207 69 

thus leading to the polynomial wu, 

Ug = D Gai (oa Y + aye Ye “f= inv" + se i) (7) 
The origin of time was chosen at the end of a coding cycle of a particular 


line. At the end of 2500 clock cycles, the cyclic process is repeated 
and the next set of calculations may be started at this instant. 


2.1.1 Lffect of Corners in the Circuit 


Corners were ignored in the polynomial calculation in the previous 
section. Considering their effects, we have 


ul = XY tty s and, uy = XY"? ty, (8a; 8b) 
from Sec. 4.3.4 in Ref. 1. However, at the gate g, the two inputs uw, 
and u, should be in phase. This condition implies that the generator G 


should generate the position a,, half a clock cycle after the instant as 
assumed in the previous calculation, leading to 


t=3 
ae oe 1 ata; 
U; = XY 2.4! = xX? as o4 t) (9) 
i1=0 
1=3 
L 93_ 5 
U, = xX”? G4 (10) 
7=0 


The element wu; has a —90 degree corner at the gate, and hence, 
i=3 
Us = bie > Ce Gates (11) 
7=0 


The constants a;; @ = 0 through 3, 7 = 3 through 5) are the same as 
their previous values. The sensors A and I should be read during the 
cycle ending at 225 + (8 — 23) = 273 clock cycles, after the generation 
of the first data position do, . 

The polynomial u; makes a —90 degree corner, and u; makes com- 
pensating +90 and —90 degree turns, before reaching P, yielding 


uh = XV ty, ; and, uf = X°Y*-u, (9a; 11a) 


*Tt can be seen that the top section of the storage loop need not have two 
independent element numbers 2 and 6; but such a numbering facilitates the repre- 
sentation and its boundary may be considered to lie at any period z located at Y 4? 
and Yoto2150+2, 
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leading to 


wo 


t= 


wie 


Us = x? a6 Youre: (12) 


=] 


= 


If incremental space positions are considered, ws after half a clock 
cycle is 


i=3 
Us = x” » ais ae (13) 
i=0 


The polynomial us has moved around two —90 degree corners before 
it reaches J, . Hence 


i=3 
, 2447 1724462 2500 24491 -7 
Uy = XL aes = hig? 
1=0 


t=3 
Ser Dy G2 ¥2 (14) 


An additional half cycle is now necessary to obtain incremental space 
positions, which results in 


i=3 
Ly7h 5001 ont 
te = X*V*us = XPM YS ae Vi (15) 
7=0 


indicating that the bubble position a), generated by G is at lp , 25003 
clock cycles from the prechosen origin of time. In this case it is at l, 
25003 cycles after its generation (having lost 4 cycle at the two —90 
degree corners of the loop), and it is going in a direction opposite to 
the direction at the prechosen origin of time, i.e., at the instant of its 
generation. 


III. APPLICATION OF THE ALGEBRA TO (7, 4) HAMMING CODE, MAGNETIC 
DOMAIN ENCODER 


3.1 Princtple of Operation 


A configuration of the encoder’ is shown in Fig. 2. The four incoming 
data bits of every code word are generated by G uniformly during 
28 ¢ seconds, where ¢ denotes the interval for one rotation of the main 
magnetic field (see Ref. 2). The data bits are accumulated in adjoining 
T-bar periods in loop 1 and gated out by g; to a duplicator D. One of 
the two resulting streams is allowed to circulate in loop 3 and the other 
is divided’’* by the generator function 


COS ae ea. 
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LOOP 3 LOOP 5 










4 OUTPUTS 
ONE EVERY 


4 INPUTS 4t (DATA) 


ONE EVERY 


3 OUTPUTS 
ONE EVERY 
4t (PARITY) 


12 Be iIVE—OR 
CIRCLED NUMBERS INDICATE ee OneG i 
ELEMENT NUMBERS. LOOP 2 


Fig. 2—(7, 4) Hamming encoder with magnetic domains and T-bars. 


There are four steps in the division cycle, each step being accomplished 
as the data stream passes through the gate g. After 4 steps and 3 circula- 
tions of the stream in loop 2 (6 periods), the remainder (the three 
parity bits) are gated by g, from loop 2 to loop 4. Meanwhile, the data 
in the loop 3 (7 periods), also having completed three circulations, is 
gated by g; to loop 5. Loops 4 and 5 have 3 periods each, and the data 
or parity bit 1s read by 8, or S, as the leading bit of the circulating 
stream in loops 5 or 4 respectively. The gates g; and g, operate iden- 
tically by diverting the bit just read by S; or S, into the annihilator A. 
Therefore, the outgoing information (7 bits) is read every 4¢ while the 
received information (4 bits) is generated every 7¢. 


3.2 The Algebra of the Circuit 


It is necessary to choose an origin of time and proceed in the time 
domain step by step as the algebra progresses. The exponent of X 
indicates the number of clock cycles between a prechosen origin of 
time and the instant of consideration. The binary values of the bit 
positions are indicated by the values of a, and the location in the circuit 
is indicated by the subscript of Y (the element number) and the ex- 
ponent of Y (the particular period in that element). 
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This circuit operates on a repetitive basis. The information is being 
continuously received at G, and the coded words are being continuously 
read at S, and 8, . The origin of time can thus be chosen at the end of 
a cycle during which the first binary position of any particular data 
block is being generated by G. 

To facilitate the representation, it is advantageous to divide the 
algebra of the circuit into a series of operations* corresponding to the 
subfunctions in the circuit. 

If the four positions of the data string generated are represented! 
and considered 21 clock cycles from the prechosen origin of time 


UL = Gains fis + Qi Y>- + Ge, Va + G31 Y;) 


} 


t=3 
ae BGs Dy aaVie. (16) 


Operation 1: Transportation of ui from generator (element 1) to loop 1 
(element 2): With their boundary located at Y? and Y$ , we have 
(from Sec. 4.3.2 of Ref. 1) 


i=3 
3 7(3-1) +3 
ui! = >» Gas ‘Ul = xe Qin Yi' ae (16a) 
i=0 
i=3 
7(3-7) 
us = Xo ai2Yo” * (17) 
i=0 


since Y2** = Y%. 
Operation 2: Looping of us in element 2 with six periods in the loop: 


i=3 7=3 
3-74 a6 -7t 
Us = xt ya i YE t) mo = xo > irs ead (18) 
~=0 7=0 


Operation 3: Gating the stream uz out of loop (element 2) to the path 
(element 3) between loop and duplicator D: If the gate g, is at Y3 , then 


i=3 
us = X*Y"-u, = X”* do QY3”, (19) 


since Y3** = Y$ when a = O while the gate g, is operational. 


DESIGN PARAMETER 1: Operating the gate g, : This gate should be 
operational for the 25th, 26th, 27th and 28th cycles from the origin 
of time. 


* This type of distinct numbering of operations is suggested when the circuit 
accomplishes complex functions. 

+ The prime (s) indicates that the polynomial as such does not represent a bubble 
stream. But after certain ensuing algebraic operations they will represent observable 
streams. 
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Operation 4: Transportation of us from gate to duplicator and its duplica- 
tion: Let the duplicator D be located at Y}. Then 


t=3 t=3 
ay. 32 (3-7), 32 (3-7) . 
Us = XxX > Nave a : } U = xX > A4¥s , . (20; 21) 
t=0 t=0 


Operation 5: Transportation of us and us to loops 2 and 8 respectively: 
If the encoder design has 9 periods in: the paths between D and loops 
2 and 3 then 


i=3 
ip = OV a = X* YS avo, (22) 
1=0 
and 
i=3 
ty = X°Y?-ug = X" YVanY?”. (23) 
1=0 


The binary positions a;; (¢ = 0 through 3 and 7 = 1, 2, 3, 4, 5, 8 and 9) 
have not undergone any transition points in the circuit. However, 
in the polynomial u;, the binary values of bubble positions a,; will 
undergo definite changes as u; is divided by the generator function. 
The effect of these changes may be represented in the algebra by dis- 
criminate use of superscripts for a;; which leads to 


Us = X"' (aos Vv: a Biek i ee fe 35 Ys). (22) 


Operation 6: Generation of the divisor polynomial Ug) (gO for the first 
tume) by the generator G, : It is seen that the instant of start of the 
generating cycle for the first bit position is not known from the pre- 
chosen origin of time. For this reason we may assume that this interval 
of time is gO and determine its value as the algebra progresses. 

Four binary positions* are generated by G, . Three cycles after the 
generation of the first position the bubble stream may be written as: 


Ugo = X” (og Ye + Gig¥e + Go, ¥e + 3, Y,). (24) 


Operation 7: Transportation of ug to gate g: If the designer has al- 
located 4 periods between the generator G, and the gate g, then’ (Sec. 
4.3.2 of Ref. 1) 


i=3 
tS ee ae, (24a) 
i=0 


* See Ref. 2 for the details of magnetic domain encoding and decoding. 

t The prime indicates that the polynomial uw, in this equation has already under- 
gone the first. step of the division cycle, and the binary values are no longer the same 
as in the previous polynomial. 
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Operation 8: Transportation of u; to gate g: If the designer has permitted 
2 periods between the entry of the loop 2 and the gate g then 


i=3 
uf = X°Y’ us; = X* > ats YE. (25) 
7=0 


DESIGN PARAMETER 2: The instant of generation of Ugo : If the gate g 
is located at Y3 and Y?, and if aos in (25) interacts with ao, in (24a), 
then they should pass through g during the same clock cycle, and the 
exponents of X associated with ao; in (25) and ao, in (24a) when the 
corresponding exponents of Y; and Y, are 2 and 4, may be equated 
yielding 48 — 3 = g0 +7 — 30r 


g0 = 36 clock cycles. (26) 


The implication of this equation is that the generator G, must generate 
Ugo (if dos is one, Ref. 2) with its bubble position a, , exactly during 
the 36th clock cycle from the prechosen origin of time. This is depicted 
in Fig. 3a. 


Operation 9: Exclusive-or operation of us in (22) and Ugo in (24): At this 
stage of the calculation it is necessary to assign binary values for the 
Qos » G15 , A>; and ag, (which are the same as data bits ao, , Gy; , @i2 and 
3 respectively) and let these be chosen as 1111 respectively. The values 
Of Aog » Gig , Mo, and az, are 1101 respectively, if the generator function 
(see Ref. 3 or 4) of the code has a form as defined in Section 3.1. 

Now us may be written as* 


i=3 
Wig Ae Do GY es (27) 
7=0 


where a}, = a°; @ a,, thus yielding a}; = 0, ai; = 0, a}, = landa;, = 0, 
SINCE Ag , Aig » A2g aNd 3, correspond to 1101. 

It can be seen that aj, is always zero and it can be dropped from 
the notation, since it is never again activated in the circuit. This leads to 


1=3 
Us1 = x » re vane (27) 
t=1 


Operation 10: Generation of the generator polynomial u, (gl for the 
second time) by the generator G, : Four binary positions are generated. 
Three cycles after the generation of the first position the binary stream 


* The second subscript 1, of w, indicates that it has gone through the exclusive-or 
operation once. 
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i=O j=1 | i=2 
IN EQUATION (17) 


















Cea AND GATING 
Y INDICATES THE OPERATION OF F DATA BITS 
Wi CHANNELLING GATES Qj, Go, 93, G4, OR Gs 

\) INDICATES THE OPERATION READ 
WS OF EX—OR GATE 


READING AND GATING 
OF PARITY BITS 





GATE 





ag a, a 





a, as] a 


t=70 80 90 


Fig. 3—The timing diagram for the various functions in the (7, 4) Hamming 
encoder. (a) The generation of data bits ao1, @11 , @21 , G31 and subsequent looping, 
gating, duplicating functions. wgo indicates the generation of the general polynomial 
by Gz. (b) The generation of wei by Gg . (c) The generation of wpg2 by Gz . (d) The 
generation of wg3 by Gz. (e) The polynomial Uy and its gating to element 10. (f) The 
parity bits [polynomial ws, , eq. (37)] and its gating. (g) The encoded data a» through 
a, and its reading by Sg and Sp 


may be represented as 
i=3 
Ui ae elts oo Gy, (28) 
t=0 


where gl is another design parameter to be determined as the algebra 
progresses. 
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Operation 11: Transportation of u,, to gate g: 
t=4 : 
tn = XT Yo aYs ”. (28a) 
i=0 


Operation 12: Looping of us, (26) once in loop 2: Loop 2 has six periods 
and after 7 clock cycles we have from Sec. 4.3.3 of Ref. 1 


t=3 
i = OY aS XI ng SX aE. (29) 
t=1 


DESIGN PARAMETER 38: The instant of generation of u,; : If the encoder 
is to function properly (see Ref. 2) then ao, in (28a) should interact 
with a;, . A calculation similar to that in operation 8 yields gl = 43. 
It is to be noted that G, generates its first binary position dp, in (28a), 
during the 48rd clock cycle from the prechosen origin of time, which 
corresponds to the generation of a), by G. (See Fig. 3b.) 


Operation 13: Exclusive-or function between ut, and u, : For the correct 
functioning of the encoder, the generator G, generates the sequence 
of four bubble positions 1101, only if the leading bubble position in 
polynomial uf, is one. In this case it can be seen that ai, is zero, so 
the generator G, generates a sequence of 4 bubble positions whose 
binary values are zero thus leading to 


1=3 

0 6-7 5 2 

Us. = x° > aa ? = Gigk Ys ) (30) 
t=1 

where ai; = @}5 ® au-y, thereby yielding aj, = 0, a3, = 1, a3, = 0, 

G5 = O Since Aog , Aig , A2z and a3, are 1101. Once again, a}, being always 

zero, can be dropped from the equation leading to 


t=4 


eg =X? a VO, (30a) 


i=2 
The last term in (30) and (30a) is the binary position a;, of (28a) with 
the exponent of X being 50 since gl was calculated as 48 clock cycles. 


Operation 14: Generation of generator polynomial u,. (¢2 for the third time) 
by the generator G, : Four binary positions are generated. Three cycles: 
after the generation of the first binary position, the binary position 
is written as 


1=3 
iS xe Ss ai.Vo° (31) 
1=0 


where g2 will be evaluated as a design parameter. 


BUBBLE CIRCUITS 1571 

Operation 15: Transportation of Ugs to gate: 
Ug = XP? > Ga: (31a) 
Operation 16: Looping Usz in (80) tn loop 2 once for 7 clock cycles: 
Ui ee eV iy SX cE a Gea (32) 


DESIGN PARAMETER 4: The instant of generation of Ug : Equating the 
exponents of X associated with the interacting terms a), in (31a) and 
@25 in (82) when they pass through the gate g we have g2 = 50 clock 
cycles. 


Operation 17: Hxclusive-or function: 


i=5 


U33 = Da >» Gio (33) 


1=3 
where a3; , Qj; , @;; may be evaluated as 101 respectively, since a,, , 
Q., and az, are 101. 


Operation 18: Generation of uss : Three cycles after the generation 
of do, we have 


ee Gale 3 ee eae (34) 
Operation 19: Transportation of Uz 3 : 

tee xe > ns (35) 
Operation 20: Looping of us3 in loop 2 

Ue = x 5 re Sama (36) 


DESIGN PARAMETER 5: The instant of generation of u,z3 : The value of g3 
can be calculated as 57 clock cycles (by equating the exponents of X 
in (85) and (86); see also Fig. 3d). 


Operation 21: Hxclusive-or function between ul, and ul, : 
i=6 
64 4 (8-7) 
Usa = xX Pe ais Ys : ) (37) 
t=4 


where ai; , @5 , @5; correspond to 111 respectively (since a;, , da, and 
a3, are 101, thus leading to the parity bits for the (7, 4) Hamming 
code with the four data bits as 111 and 1. 
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Operation 22: Gating of parity bits from loop 2: Now the parity bits 
at Y?, Y2 and Y? can be channeled out of the loop 2 into the element 6 
by the action of the gate g, located at Y;. This function may be repre- 
sented as 


Il 


a=6 
OV aga x Dae 
i=4 
= X(aysV5 + ds5V¥s + a5sV3), (38) 
since Y2** = Y% when a > 0 thereby leading to 
= KG e os ight + Age Ya: (39) 


DESIGN PARAMETER 6: The operation of the gate g. : The gate is located 
at Y? and Yé and it can be seen that a4, , @5, and ag, reach Y§ when the 
exponent of X is 65, 66 and 67, indicating that the gate must operate 
during these three cycles to channel the bits for loop 2. 


Operation 23: Looping of uy three times:* The data stream uy has been 
circulating in loop 3 for (8 X 7) clock cycles leading to 


t=3 t=3 
-i 62 3-7 
>= > ee >> dave t) mod 7 = XxX > aiYVe . (40) 
7=0 t=0 


Operation 24: Transportation of uy to gate gs : This gate is located at 
Y3 to match the location of the gate g. at Y} (see operation 22). If 
the gate operates at the appropriate time for four clock cycles, then the 
polynomial representing the stream in the path between loop 3 and 
loop 5 is 


t=3 
6 a XY? uy as x os OY er, (41) 


t=0 


But Y3** = Y,¢ during gating of g; and therefore we have 
Uy = Dus x a; yo" ae 


DESIGN PARAMETER 7: The wherdlian of gate g,; : The polynomial wo 
in (41) indicates that the gate g,; should operate when the bubble 
positions a) , @, , @ , and a; are at Y{, . Further, when the exponents 
of X are exactly 64, 65, 66 and 67, the four binary positions ay , a , 
and a3 are in the gate’g; . Hence, the operation of this gate must coincide 
with the 64th, 65th, 66th and 67th clock cycles from the prechosen 
origin of time: 

* This operation takes place during the 26 clock cycles allocated for operations 6 
through 22 for the polynomials ws, ws1, Ws2, 53 and wo. 
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The single clock cycle between the operation of gates g; and gs 
(Fig. 8e and f) is due to the incomplete fourth rotation of the parity 
bits us . These are gated out by g; after 3 rotations and after 4 steps 
of the division. (See Ref. 2.) 


Operation 25: Transportation of data and parity bits to loops 5 and 4: 
If there are 4 periods in elements 10 and 6, then 


i=3 
Un = X*Y*-uU19 = x” > G:Var 5 (42) 
i=0 
and 
i=6 
uk = X*Y*-ug = X" Da, YP. (43) 
i=4 


DESIGN PARAMETER 8: Reading of data bit a) and its gating by gs into 
the annthilator A: If the sensor Sy, is located at Yj, , then ay will be at 
Yi, at X°, indicating that the sensor should be read during the 69th 
clock cycle from the prechosen origin of time. If the gate g; is located 
at Y?, , then ay is at Y?, at X’° and it should operate during the 70th 
cycle. After the gating of ad) by g; , we have 


t=3 


Un = Ge ae ie (44) 


i=1 


Operation 26: Looping for 2 clock cycles: After 2 clock cycles, a, will 
be at Yj, , and it can be read again, since 


1=3 
73 (5-121) mod 3 
lig ae SY (45) 
t=1 


DESIGN PARAMETERS 9, 10, 11: Reading and gating of a, , a2 and az 
by Sq and gs respectively: The sensor 8, is read during the 73rd clock 
cycle and gate g, operates during the 74th clock cycle from the origin 
of time. More operations of the type 27 indicate that Sy should read a, 
during the 77th clock cycle, and g; should gate a, during the 78th 
clock cycle, Sg should read a; during the 81st clock cycle, and g, should 
gate a; during the 82nd clock cycle. These details are plotted in Fig. 3g. 


Operation 27: Looping! of parity bits in loop 4: After 10 clock cycles 
the parity bits in loop 4 may be represented as 


* It is no longer necessary to carry a second subscript for a, since none of the 
binary values in any of the polynomials change in the remainder of the circuit. 

t This operation takes place during the 10 (i.e., —2 for ao, and 4 for a1, a2 and a3 
each) clock cycles allocated for reading, gating and looping of the four data bits 
Qo, G1, Ae and a3. 
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0 d 3 d3 0 d 3 
Ur = xo Gyo? ae GY + As ; mo ) 


= aa Pe oe + ie ee ++ GaN a). (46) 
Operation 28: Looping of u; : After 4 clock cycles* 
Uz; = X”*(asY1 + asY? + a,Y3). (47) 


DESIGN PARAMETERS 12, 13 AND 14: Reading and gating of a, , a; and 
a, : It is seen that a, is at Y; (location of 8,) during the 85th clock 
cycle. The gate g, at Y7 should divert a, during the 86th clock cycle. 
More operations of type 28 indicate that S$, should read a; during the 
89th clock cycle, and gs should gate a; during the 90th clock cycle, 
S, should read a, during 93rd clock cycle, and g, should gate a, during 
the 94th clock cycle. These details are also plotted in Fig. 3g. 

A block of data is thus completely processed by the circuit after 
the 94th clock cycle, and the position of every bit of information may 
be accurately predicted during any prechosen clock cycle. The design 
parameters are also accurately determined by the analysis. 


VI. APPLICATION OF THE ALGEBRA TO 16 LINE LDC MAGNETIC DOMAIN 
LINE SCANNER 


4.1 Principle of Operation 


The circuit for the line scanner is shown in lig. 4. Sixteen inputs, 
1 through 16, carry telephone line currents in a loop for each line. 
When the line current is not sensed (i.e., on-hook status) during a scan 
interval, a bubble in the position (1-2), (8-4), --- (81-32) of a 75- 
position loop is moved from the outer periphery to the inner periphery 
of the loop, thus effecting a Lateral Displacement Coding. 

The laterally displaced bubble positions are moved under two sensors 
S, and S, . The spacing between these is arranged to sense the status 
of a particular line at two instants of time. When the sensor 8, is sensing 
the contents of a certain cell coded by the circuit in line 7, (j = 1 through 
16) at an instant ¢, then the sensor 8, is sensing the contents of the 
cell Y coded by the current in the same line j at ¢ + 6. (The value of 6 
is the duration required to move the bubbles from sensor 8, to 8, and 
is also the scanning interval of the lines.) When there is no discrepancy 
between the readings of S, and §, , then there is no change in status 
of the line (on hook or off hook) during the scanning interval and vice 
versa. The status of the jth line is coded during the movement of the 


* It is interesting to note that the extra clock cycle between the gating of g3 and ge 
(design parameter 7) is really necessary for correct functioning at this stage. 
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Fig. 4—The schematic diagram of 16 Telephone-Line-Scanner. 


bubble across the cell from a position (27 — 1) to 27. A conductor 
carrying a single phase current around the loop propagates the bubble 
positions around its periphery at a uniform speed. Further, when the 
bubble positions pass through the cell 74-0, they are all positioned 
towards the outer periphery of the loop, and are reset ready to be coded 
again. 


4.2 Algebra of the Circutt 


The origin of time may be chosen at the start of a coding cycle which 
repeats every 36 clock cycles. The lateral displacement coding occurs 
during the single clock cycle immediately after all the cells 0 through 35 
contain bubbles at the outer periphery of the loop. The start of this 
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cycle would then constitute the origin of time, and the polynomial 
representing the stream in this section of the loop is* 


1=3 


Uoo = x aye (48) 


7=0 


on 


where 7 = 0 denotes the leading bubble position in the stream; and a; 
is the uncoded status of the bubble position. Next consider the bubble 
stream in the section 36 to 71 of the loop. This stream carries the status 
of each line from the last cycle, and the polynomial representing this 
stream at the prechosen origin of time is 


1=35 
tho = X° a. ¥"™, (49) 
1=0 
where 2 indicates the bubble position in the stream. Similarly, the 
polynomial wu, between 72 and 74 may be written as: 


1=35 


ig Sa > Gee (50) 
t=33 
The value of 7 ranges from 33 to 35 since there are only three locations 
72, 73, and 74 to accommodate the last three bubble positions of the 
stream coded during their own coding cycle. The contents of the entire 
loop at the origin of time may be written as: 


Uo = Uoo + U1o + Uso 


1=35 1=35 1=35 
xy ay * +4 SS Ci fae + > ae) (51) 


t=0 t=0 t=33 


Now consider the movement of the bubble stream represented by (51); 
for m (m S 36) clock cycles, the resulting polynomial according to 
Sec. 4.3.2 of Ref. 1 is: 


t=m-1 1=35 


Ul = xr( ys acy ate > oe Games 
0 7=0 


i= 


1=35 1=35 
as > Cr? + Pe wee) 


=0 t=33+m 
or 
t=m—-—3 


= Uor + Ui + Ua + Uar , (52) 


where Uo, represents the first m binary positions of wo. being generated 
in the section 0 through 35 of the loop; u,, results from the translatory 


on 


* There is no need to subscript Y since there is only one element (the loop) in 
the circuit. 
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movement (Sec. 4.3.2 of Ref. 1) of woo in (51), and ue; results from the 
movement of wo in (51). The lower limit of 7 in u., should be chosen 
according to the value of m. When m is less than 3, the value of 2 equals 
0 is appropriate, since the leading bubble position a») is within the 
maximum number of binary bubble positions (i.e., 74) in the loop. 
When m exceeds three, the leading bubble ay of uw; in (51) is transformed 
as the fourth term of wo, in (52). The polynomial ws, is the transforma- 
tion of wea when m is less than 3. When m exceeds two, ws; drops out 
of (52), since the lower limit of 7 exceeds the upper limit of 35. 

Examine the polynomial u, in (52), when m = 0, wo; drops out of 
the equation, and w,, , us, and u3, will assume the roles of the poly- 
nomials Ugo , Uio and Uo in (51). Next observe the polynomial uw, in 
(52); when m reaches a value of 36, uor , Ui1 and U2, Will assume the 
roles of Woo » Uo and wo in (51), and uz, drops out of the equation. In 
essence, we have two cyclic processes taking place simultaneously, 
the first one being in the time dimension and repeating every 36 clock 
cycles, the second one being in the spatial dimension and repeating 
every 75 periods. The effect of the first cyclic process may be eliminated 
by always considering m as (m mod 36). The effect of the second one 
may be eliminated by always considering the exponent (e) for Y as 
(e mod 75). 

Table I relates the values of the exponents of Y for various values 
of 2 and m in the polynomials wo; , Wi1 , Wer , aNd uz, . It also indicates 
the locations of the first and last bubble positions of streams repre- 
sented by these polynomials. 

4.4 Implication and Use of the Representation for the Line Scanner 


4.4.1 Prediction of Bubble Positions 


Consider the tenth (¢ = 10) bubble position in a data stream coded 
for twenty (m, = 20) cycles and propagated for ninety (m = 90) cycles. 

The initial position under consideration is a,).X°’Y'. Table I indicates 
that this term exists in u,,; with a value of 1 = 45 yielding the location 
of this position. When this location is propagated for 90 cycles, the 
new position is Y“°t°? ™-475 ie. Y°°: and the corresponding value 
of m* is (20 + 90) mod 36, i.e., 2 eyeles. The bubble position is then 
a,X°Y. The only positive value of 7 which satisfies the constraints 
on the exponents of both X and Y is 13; and the individual term de- 
noting this bubble position lies in the polynomial uw, in (52). This 
implies that if the original position in uw, is @iX°°Y*’, then after 90 


* When the exponent of Y is less than (m — 1), it should be concluded that the 
term is in wo1, (see Table I) and is in the dead interval of the circuit. 


TaBLE I—LocaTION oF BUBBLE STREAMS 


Limits of the exponents of Y in the polynomials uy, U4, Ux, and Us, in (52) 














m1 85 35 35 
m or Uo = xm > a;Yyr-\-i 1 = Nm > a, Yttm—7t Uy = Xm > a,;Yitm-i Usa = xm es a,;Y'or-t 
m mod 36 0 0 0 or m—3 83+m 
0 = 35 to 0 71 to 36 
(2 = Otoi = 35) (i = 0 toi = 35) 74, 73, 72 
1* 1=0 36 to 1 72 to 37 74, 73 
(i = 0 toi = 33) (i = 0 toi = 33) 
2 1,0 37 to 2 73 to 38 74 
3 2 to 0 38 to 3 74 to 39 — 
4 3 to 0 39 to 4 74 to 40 — 
(@ =m — 8) 
20 19 to 0 55 to 20 74 to 56 — 
m (m — 1) to 0 35 + m tom 74 to 36 + m — 
(t =m — 8) 
34 33 to 0 69 to 34 74 to 70 —= 
35 34 to 0 70 to 35 74 to 71 — 
36+ ee 35 to 0 71 to 36 74, 73, 72 








* 0-1 constitutes the coding cycle 
} See m = 0 values 
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clock cycles the final position in wu; is a;,;X°Y” serving as the thirteenth 
bubble position in the data stream. I'urther, during the coding cycle 
i.e., as the exponent of X is changing from 0 to 1, the location of the 
initial position is at Y*°-*°, ie., at Y°° serving as an active bubble 
carrying the status of the 13th line in Fig. 4. During its next coding 
cycle, the same bubble position is at Y°°-’, or at Y°* corresponding 
y**-** or at Y*’ during coding, serving as an inactive bubble position 
between lines 1] and 12 in Fig. 4. 


4.4.2 Detection of the Nonexistence of Bubbles 


For the correct functioning of the line scanner, all the bubble positions 
should carry bubbles. Sixteen lines are actively used; and the status 
of these lines is carried by bubbles in positions 1, 3, --- 31. The bubble 
positions 33 and 35 always carry the status of two fictitious lines (on 
on-hook and next off-hook) to check the correct operation of the overall 
magnetic and electronic circuitry. Generally, it is also desirable to 
check if all bubble positions do carry bubbles by sensors 8, and 8, 
which are capable of detecting only the off-hook status of lines. Such 
an inspection can be effected when each of the bubble positions is 
arranged to periodically occupy the position 35, which should always 
carry the off-hook status. All the bubbles are moved to this status as 
they traverse the position 74. 

Examine the bubble position at Y°*’ just prior to coding. After 36 
cycles (i.e., next coding) the bubble position now at YS °° ™?” ie. 
Y’* will occupy Y°*’. In general, after n coding cycles, the present 
position Y°?-**™ ™475 will occupy Y°*. The exponent of Y generates 
a series 35, 74; 38,2; --- , (85 + 3(n)/2), 74 + (8(m — 1)/2) mod 75; 
-++ , 32, 71; repeating every time n reaches 25. This indicates that the 
bubble positions now occupying Y*, Y° --- , ete., Y", Y’ --- , ete, 
never occupy Y°’ at any finite value of n. To eliminate this condition, 
the electronic circuitry may be programmed* to delay the coding by 
one clock cycle every 25 coding cycles. This leads to a new location 
series: 35, 74, 82, 2, --- , 82, 71; 36, 0, 39, 38, --- , 83, 72; 37, 1, 40, 
4, --- , 34, 73; 38, 2, --- , ete. Alternatively, if the loop is designed 
with 73, 109, 145, ete., periods, then the need for building additional 
delay circuits will not be necessary. With 73 periods, every bubble 
position will be located at Y*’ every 73 coding cycles, or every 2628 
(i.e., 73 X 36) cycles and so on. 


* The general concept of shuffling periodically was suggested by D. Denburg. 
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V. CONCLUSIONS 


The polynomial algebra is a flexible mathematical tool available 
for the step by step design of conceived circuits, and for the sequential 
verification of their operation. Various design parameters may be 
calculated accurately. 

When the operations of numerous circuits are to be synchronized, 
the algebra provides an excellent insight into their combined functioning. 
The effect of errors or defects of certain sections of the overall circuitry 
may also be accurately analyzed by the algebraic modeling of the 
circuit operation. 
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XYTOLR-A Computer Program for 
Integrated Circuit Mask Design Checkout 


By MICHAEL YAMIN 
(Manuscript received March 7, 1972) 


XYTOLR is a computer program which checks integrated circuit mask 
designs for compliance with design rules such as those which require 
minimum clearances between diffused and metallized areas. The input 
ts the same machine readable mask description that is used for automatic 
mask art work generation. Clearance violations, if present, are reported 
in graphic and printed form with pointers to the erroneous statements 
in the user’s input. 


I. INTRODUCTION 


The increasing size and complexity of large-scale integrated circuit 
mask designs has created a demand for computer aids to the circuit 
designer. Without such aids, the design of many circuits would be not 
merely difficult, but practically impossible. Perhaps the best established 
design aids in this field are computer art work generation systems, 
which accept as input numerical descriptions of mask layouts, and drive 
high-precision graphical output devices which create the mask masters. 
XYMASK,*’’ used at Bell Laboratories, is such a system. 

The use of a computer art work generation system implies that there 
exists a complete and unambiguous machine-readable description of 
the masks required to produce an integrated circuit. In principle, this 
information, together with the physical properties of the substrate 
and diffused layers, is sufficient to predict any desired property of the 
completed circuit. The availability of computer programs to predict 
and analyze important properties of circuit designs in advance of their 
fabrication would be a major aid to the circuit designer. 

One early phase of the design process for which a computer aid is 
highly desirable is that of checking a mask set, before it enters pro- 
duction, for compliance with design rules such as required clearances 
between areas delineated by the same or different masks. Deviations 
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from these rules may result in loss of manufacturing yield or early 
device failure. Visual checking of mask art work for such tolerances 
is not only tedious and time-consuming, but has frequently been found to 
miss one or two significant errors in the mask design. The exhaustive 
patience of the computer is better suited for such important but 
mechanical verification procedures. This paper describes a computer 
program called XYTOLR, which tests mask designs automatically 
for tolerance errors. The data structures and other components of 
XYTOLR can be used to construct other mask analysis programs. 


II. USER INTERFACES-—INPUT AND OUTPUT 


Computer programs, in general, should require the minimum neces- 
sary input from the user, and should accept it in a form which is simple 
and natural for the user to write. The output should include all the 
information which the user might desire, but indexed and organized 
so that specific items can be located without extensive searching. An 
attempt has been made to follow these principles in the design of 
XYTOLR. 


2.1 Mask Design Input 


The primary input to the tolerance-checking program is the geo- 
metrical description of the mask set in XYMASK'” input language. 
This is the same input deck (except for a few command statements) 
that 1s used for mask art work generation via XYMASK. Since such 
decks can run to several thousand cards, it is unreasonable to ask the 
user to make any general alterations for the benefit of XYTOLR. 
To preserve language compatibility, a version of the XYMASK language 
processor is incorporated in XYTOLR to read this input. 

In the XYMASK language, the mask designer defines closed plane 
figures composed of line and circular are segments. Each such defined 
figure is called an “atomic,” and its definition includes a symbolic 
name (mask level name) which assigns the atomic to a specific mask 
of the set being designed. An atomic may also be given a label (atomic 
name). By the use of statements referring to this label, an atomic 
may be translated, rotated, reflected, and repeated many times at 
different locations in the mask. Groups of atomies may be combined 
into compound structures called “clumps,” which may also be labelled 
and manipulated like atomics. A ‘‘clump” may include atomics assigned 
to different mask levels. In this way, devices such as transistors which 
require shapes on several masks can be designed and manipulated as 
units. 
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While XYMASIX permits atomics to be defined without labels 
(atomic names), XYTOLR requires a name for each atomic. A name 
is therefore generated for any atomic not named by the user. Such 
names are of the form ATnnnnn, where nnnnn is a number pointing 
to the location (line number) of the atomic definition statement in the 
XYMASK deck listing. The user thus can find the statement to which 
the generated name refers. 

The structure described by a XYMASIK< deck consists essentially 
of many geometrical figures distributed among a number of mask 
levels. Each figure is called an ‘‘occurrence.”?’ XYTOLR assigns a serial 
number to each occurrence. A single atomic definition may result in 
many occurrences, as a result of repetitions. 


2.2 Tolerance Test Input 


XYTOLR, like XYMASK, is technology-independent, dealing only 
with the geometry of the mask design, and knowing nothing about 
design rules or the device functions of the various masks. The user 
must provide input which describes the ‘tolerance tests’? he wishes 
the program to conduct. A “tolerance test’? asks for a report of all 
instances in which the clearance between elements on specified mask 
levels is less than a specified quantity. The operation involved in a 
tolerance test is as follows: enlarge each occurrence on a designated 
mask level by some quantity and report all new connections which 
develop, and old connections which are broken, between these occur- 
rences and the occurrences on another designated mask level. The 
enlargement quantity is generally a little less than the specified mini- 
mum clearance, so that clearances just at the minimum are not reported 
as errors. If checking of clearances within a mask level is desired, the 
operation is: enlarge each occurrence by a little less than half the desired 
clearance and report all new connections within the mask level. 

A simple, free-form language is used to specify tolerance tests. For 
example, suppose one wishes to assure that every emitter in a bipolar 
integrated circuit lies no closer than 0.5 units to the boundary of the 
base diffusion which surrounds it. (The units are the same as those 
used in the XYMASK< design deck.) The statement 


EMITDIF +.5 BASEDIT 


results in every emitter being enlarged by a little less than 0.5 units, 
and any new contact with the base diffusion boundaries being reported. 

Suppose one wishes to check that every occurrence on a metallization 
mask clears every other occurrence by 0.3 units. The statement 
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/ THINMET +.15 


results in every occurrence on the THINMET level being expanded 
by a little less than 0.15 units, and any new connections within that 
level being reported. 

A negative number results in contraction of occurrences, and may 
be used to test for minimum overlaps. 

The tolerance test language has the capability of describing quite 
complicated test patterns within a statement. It includes control words 
which enable the user to modify the normal handling of a tolerance 
test by XYTOLR: for instance, expand by exactly the quantity re- 
quested rather than a little less. Any number of such statements may 
be submitted, and each will induce a tolerance test and generate a 
report. 

In practice, a mask designer working on a particular type of circuit 
will have developed a standard test protocol to be applied to all of 
his designs. Thus, specifying a tolerance test sequence will consist 
simply of selecting the appropriate tolerance test input deck for the 
circuit type being tested. 


2.3 Masklevel Connectivity Input 


The user may know that tolerance tests between certain mask levels 
will never be made, either because these mask levels can never be in 
physical contact (for example, a metallization and a buried collector 
in bipolar technology) or because there is no current interest in such 
a test. XYTOLR can accept input that tells it which mask levels are 
to be considered capable of connection, and which are not. It is never 
necessary, in principle, to supply such input, but, for large mask designs, 
considerable savings in running time and memory occupancy can result 
if it is supplied in advance of all tolerance test input. 


2.4 Output 


Both printed and graphical reports are generated by XYTOLR. 
For each tolerance test, a plot is produced of all occurrences involved 
in violations of the tolerance rule. Each plotted occurrence is identified 
with its serial number. Comparison of this plot with drawings of the 
individual masks (provided on the same scale) facilitates locating the 
errors in the mask design. 

The drawing of tolerance errors is accompanied by printed output 
identifying each pair of occurrences involved in a violation. The occur- 
rences are identified by atomic name, mask level name, occurrence 
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serial number, and coordinates. This information makes it possible 
to find the offending statements in the XYMASIX input deck. 

Self-intersecting shapes, which are not permitted by some precision 
plotters although they are accepted by XYMASK, are also reported. 

The program has been designed to err on the side of strictness in 
reporting tolerance violations. Marginal cases are reported as violations, 
leaving it up to the user to decide whether he considers them acceptable. 
If no violations are found, that fact is reported explicitly. 

Tigures 1 and 2 are examples of XYTOLR output. Figure 1 shows 
a simple metallization pattern. The clearance desired between metallized 
areas was 0.3 units. A number of occurrences in violation of this rule 
were deliberately inserted into the XYMASI< description of the mask. 


' | t t t 
DEMONSTRATION OF XYT@LR @CCURRENCES ON MASK LEVEL N19 





Hg. A simple metallization pattern. Two shapes have been inserted in 
violation of the required minimum clearance of 0.3 units. 
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A tolerance test was induced by the input statement 
/ N19 +0.15 NODEXX. 


NODEXX is a control word requesting that new connections 
(violations) involving occurrences already indirectly connected through 
other occurrences (i.e., members of the same node) be plotted separately 
from new connections between occurrences not so connected (i.e., 
members of different nodes). Figure 2 is the graphic output illustrating 
the latter set of violations, which are the ones of interest in this case. 


III. PROGRAM STRUCTURE 


3.1 Data Structures 


The output of the XYMASK language processor is a sequential file 
or tape containing an explicit description of each occurrence in the 
mask set. Normally, this file is used to drive a precision plotting device. 
XYTOLR, however, needs parallel access to many occurrence descrip- 
tions and therefore must have a randomly accessible data base, in core 
memory, or on random access dise or drum. Since XYTOLR has been 
implemented on a computer (HIS 6070) with a large core memory, 
its data structures are held entirely in core. The principles of operation 
would not be significantly different for virtual memories or paged 
memories, or for random access external storage. 

The information contained in the XYMASK output file (as generated 
by a somewhat modified version of the XYMASK language processor) 
is converted to a data structure, the unit of which is the “occurrence 
block.”” Each occurrence block describes one occurrence. The length 
of each occurrence block is 11 + 2 + 6A core locations (36 bit words), 
where L is the number of line segments, and A is the number of arc 
segments, composing the occurrence. All geometric information is 
stored as single precision floating point numbers. 

Occurrence blocks are stacked sequentially in memory. Each occur- 
rence block contains, along with line and arc descriptions, a header 
containing the extreme coordinates of the figure (circumscribed rec- 
tangle), its serial number, an index defining its atomic name, and 
pointers which serve to organize the data structure. By means of these 
address pointers, each of which is the address of an occurrence, the 
occurrence block is made a member of three linked lists. By following 
one list (“main list’’) from the beginning, every occurrence is encountered 
in order of increasing x (min), which is the minimum z-coordinate in 
the figure (including projecting arcs). By following another list, every 
occurrence of the same atomic name is found. The third list pointer 
is the “intersection ring pointer,” by which all occurrences which 
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DEM@NSTRATION OF XYTOLR TOLERANCE TEST 1 NOTIN NODE 


Fig. 2—Output of an XYTOLR test on the pattern of Fig. 1. The input statement 
was / N19 +0. 15 NODEXX. 


contact one other directly or indirectly are ultimately chained into a 
ring, or circular list, which represents a node. 

The data structure also contains, among other things, a table of 
atomic names and of mask level names. Cross-referencing among these 
tables and the occurrence blocks facilitates finding the address of a 
desired occurrence block or the properties of an occurrence the block 
address of which is known. 

One location within the data structure contains an important number 
known as the “critical distance.”” Two points separated by less than 
this distance are considered coincident. The critical distance, while 
easily variable, is generally set as 0.0001 times the overall mask size. 
The quantity by which tolerance test enlargement parameters are 
modified, as described in Section 2.2, is 1.1 times the critical distance. 
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3.2 Subroutines 


The working tools required to perform analyses of the data structure 
are subroutines which generate, manipulate, and display occurrences, 
given their block addresses. Using these subroutines, along with the 
crossreferencing provided by the various pointers and tables, a pro- 
grammer may construct systems which operate on mask levels or on 
the entire mask set without concerning himself about the detailed 
structure of occurrence blocks. Manv such subroutines are contained 
in XYTOLR, performing both simple and complex functions. Two 
of particular importance are CONTAC and NLARGR. 

CONTAC tests a pair of occurrence blocks (call them “‘a’’ and ‘“‘b’’) 
for connectivity. It returns one of four numerical values, each cor- 
responding to one of the four possible classes of connectivity of a and b: 
disjoint, a inside b, b inside a, and overlapping (intersecting boundaries). 
CONTAC first checks the extreme coordinates of the two occurrences 
to see whether they are clearly disjoint. If not, the test for overlapping 
boundaries is made: each side of a is solved analytically with a side 
of b to see if there is an intersection within the bounds of both sides. 
(This analytical solution is returned by a subroutine called INSECT, 
which handles lines and ares interchangeably.) If no such intersection 
is found, the ‘‘inside” test is performed: a line segment is generated 
from a point on the smaller occurrence (in area) to a point outside the 
larger one, and the number of intersections of this line segment with 
the larger occurrence boundary is found. If this number is odd, the 
smaller occurrence must be within the larger one; if even, they are 
disjoint. 

NLARGR expands or contracts occurrences. Given the address of 
an occurrence block and a positive or negative expansion distance, 
it generates a new occurrence block representing the expanded occur- 
rence. In its normal mode of operation (‘‘square mode’’), it moves each 
occurrence side perpendicular to itself by the specified distance and 
extends it to intersect its relocated neighbors. NLARGR has a more 
precise, but more expensive, mode of operation (“‘round mode’”’) which 
inserts arcs in outside corners so that the new outline truly represents 
all points equidistant from the original. 


IV. OPERATION OF XYTOLR 


XYTOLR checks a mask design for tolerance and clearance errors 
by the systematic application of subroutines NLARGR and CONTAC 
to the data structure representing the mask design. If the connectivity 
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value returned by CONTAC for a pair of occurrences is changed when 
one or both of the occurrences is expanded by the prescribed clearance 
distance, a violation of the tolerance rule is indicated. 

XYTOLR constructs an essentially complete connectivity matrix 
for all the occurrences in the mask set. The data structure is then 
modified in accordance with the user’s instructions so that all occur- 
rences on the appropriate mask levels are enlarged or contracted; 
each matrix element is recomputed and those elements which differ 
from the original are flagged as errors. 

In principle, a complete connectivity matrix for a mask set would 
have as many elements as the square of the number of occurrences 
in the mask set—a million matrix elements for a thousand-occurrence 
mask set. The problem of storing such a large matrix is simplified by 
the following considerations: 


(i) The matrix is diagonally symmetric and the principal diagonal 
is of no interest. 

(27) Each matrix element may have only four values (the four types 
of connectivity) and therefore may be represented by a two-bit 
binary field. 

(iit) The matrix is generally sparse; that is, most of its elements 
are zero, and it may usually be reduced to a set of submatrices with 
a much smaller total number of elements. 

(tv) The manner in which the tolerance test is conducted permits 
the matrix to be scanned in a predetermined order rather than 
randomly; therefore, it may be kept in peripheral storage with only 
relevant elements in core at any time. 


XYTOLR first constructs the occurrence-block data structure from 
the modified XYMASK output file. At this stage, the main list pointer 
of each occurrence block points to the next sequential block, and the 
intersection ring pointer of each occurrence block points at the block 
itself. A sorting program rearranges the main list pointers so that 
occurrences can be followed in order of increasing x (min). 

Next, the “nodes” of the data structure are identified. A node is 
a set of occurrences each of which is connected to every other, either 
directly or indirectly through other members of the set. (Two occur- 
rences are considered ‘‘connected”’ if they are not disjoint.) If the data 
structure represented a metallization pattern, these nodes would be 
the electrical nodes; the term is here generalized to mean any set of 
interconnected occurrences. To find the nodes, each occurrence is 
tested for contact with other occurrences using CONTAC. The first 
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occurrence in the main list is tested against each subsequent occur- 
rence, similarly the second, etc. The list need be followed in each case 
only until an occurrence is found the x (min) of which is greater than the 
x (max) of the test occurrence. Geometrically, each occurrence is tested 
against all those which impinge on a vertical strip defined by its mini- 
mum and maximum z-coordinate. 

When two occurrences are found to be in contact (i.e., not disjoint), 
their intersection ring pointers are interchanged. This operation is 
sufficient to combine two isolated occurrences into a ring (circular list), 
to add an isolated occurrence to an existing ring, or to combine two 
existing rings. (Occurrences already members of the same ring are not 
further tested against one another.) When the operation is complete, 
each node is represented in the data structure as a separate circular 
list. Isolated occurrences are ‘‘one-membered”’ rings. 

Next, the intersection matrix is generated. Since an occurrence 
must be disjoint with any occurrence not a member of its node, only 
matrix elements relating members of the same node need be computed. 
Clearly, the more nodes, and the fewer occurrences per node, the fewer 
matrix elements need be computed. At this point, the nodes and their 
submatrices are stored on two sequential files, the ‘‘node file’ and the 
“matrix file.’”’ The data structure is scanned for nodes; each time one 
is found, its occurrences are written out sequentially; as each occur- 
rence is written, its connectivity values with subsequent members of 
its node are determined and written out as a sequence of two-bit fields. 
The files thus generated are used for all subsequent tolerance tests. 

Each tolerance test statement in the user’s input is decoded and 
the appropriate enlargement for each mask level determined. For each 
tolerance test, the data structure is reconstructed using occurrence 
blocks read back from the node file. As each occurrence block is read, 
its mask level is determined and the required enlargement applied, if 
necessary, using NLARGR. When the new data structure is complete, 
the connectivity values relating occurrences in each node are recomputed 
in the same order as that in which they appear on the matrix file and 
compared with corresponding values read back from that file. Any 
failure of these values to match indicates a tolerance test violation, 
and output routines are called to add descriptions of the pair of of- 
fending occurrences to the output report and write the occurrence 
blocks on a graphics file for later plotting. Finally, each occurrence 
is tested for connectivity with occurrences not members of its node; 
any contact indicates an error. 

The generation of the original noded data structure and of the node 
and matrix files, which are used in all tolerance tests, incurs a con- 
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siderable overhead expense for large mask sets. The user is therefore 
given the option to have all this information saved on a magnetic 
tape as soon as it is generated. The data structures and files can be 
reconstructed efficiently from this tape if further tolerance tests are 
necessary, or in case of abnormal termination of the computer run 
before all tolerance tests have been completed. 


V. EXTENSIONS OF THE XYTOLR SYSTEM 


The program XYTOLR is built about a data structure and a set 
of subroutines which generate, manipulate, and display elements of 
this data structure. These subroutines can be used to construct mask 
design and analysis programs other than XYTOLR. The programmer 
need never concern himself with the detailed structure of the occur- 
rence blocks or with computations involving the fundamental geometric 
data, but may interact with the data structure almost entirely through 
subroutine calls which provide the information or perform the actions 
logically required by his program. 

The program ENMASK, assembled mainly from XYTOLR com- 
ponents in this way, creates new mask level descriptions from X YMASK 
originals. The new descriptions, which are punched in XYMASK 
input language, represent mask levels on which each original shape 
has been expanded or contracted by a given quantity, so that it sur- 
rounds, or is surrounded by, the original shape. ENMASK can also, 
if desired, merge connected shapes, replacing each connected set of 
occurrences by the outer boundary of their union. This program provides 
a simple solution to certain mask design problems. It utilizes an occur- 
rence-manipulating subroutine called OVRLAP,’ which, given two 
occurrence blocks, can generate new occurrence blocks representing 
each simple closed figure resulting from their superposition. 

D. G. Schweikert has used XYTOLR components to construct a 
system called ICSORT, which can “decode” an XYMASK mask set 
description representing a bipolar integrated circuit, automatically 
identify the transistors, diodes, and resistors, and generate a node list 
topologically equivalent to the circuit diagram. The designer can thus 
assure himself that the masks he has designed in fact represent the 
desired circuit. 


VI. IMPLEMENTATION 


XYTOLR has been implemented on a Honeywell Information 
Systems (formerly General Electric) 6070 computer, under the GECOS 
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operating system. The component programs of XYTOLR have been 
written principally in FORTRAN IV. Assembly language subroutines 
have been used sparingly to perform functions unique to the HIS 6070. 
Some FORTRAN code involving character manipulation (for titles, 
etc.) has been designed on the assumption of a six-character machine 
word and would have to be modified to execute on a computer with 
a different internal character representation. Graphic output may be 
obtained on such devices as the General Dynamics SC-4060 microfilm 
plotter. Calls to the Bell Laboratories plotting subroutines have been 
centralized for easy replacement. 

XYTOLR occupies a core region the size of which it varies dy- 
namically according to the current size of the occurrence-block data 
structure. The size of the region is usually limited by system con- 
siderations to about 96000 36-bit words. Of this, the program (which 
is heavily overlaid) occupies about 20000 words, and the rest is available 
for data. Under these conditions, 2800 to 4000 four-sided figures, or 
2200 to 2800 eight-sided figures, can be accommodated, the exact number 
depending on the manner in which the mask set design is organized. 

Running time of XYTOLR depends on the size and complexity of 
the mask design and on the tolerance tests requested. Running 
XYTOLR on the simple metallization pattern of Fig. 1 took 62 seconds 
of processor time, including generation of tape to drive the microfilm 
plotter. The tolerance test itself took 18 seconds, and additional tests 
would add about this cost apiece. Extensive checking of large and 
complex mask sets has on occasion required as much as an hour or 
more of processor time. The computer charges necessary for such 
checkout are considered minor compared to the overall cost of the 
mask set. 


VII. CONCLUSIONS 


A computer program called XYTOLR tests integrated circuit mask 
designs, described in the XYMASK input language, for conformity 
with tolerance rules which require minimum clearances between dif- 
ferent regions and minimum overlaps between regions intended to 
be connected. Tests are specified in a simple input language. Graphical 
and printed output is produced which identifies all deviations from 
the specified rules and facilitates their correction. In the HIS-6070 
computer environment at Bell Laboratories, mask sets with as many 
as 3000 to 4000 individual shapes may be processed. 

Mask designers have found XYTOLR to be a valuable tool for final 
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mask checkout before production approval is granted. It is especially 
impressive when the program finds one or two demonstrable errors 
in a mask set which had been considered completely checked and 
ready for production. The program has proved reasonably easy to use 
and the output reasonably easy to interpret. 

The techniques, data structures, and subroutines which constitute 
the XYTOLR program may be used as the fundamental building 
blocks from which other mask analysis and design programs may be 
constructed. 
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Derivation of All Figures Formed by the 


Intersection of Generalized Polygons 


By MICHAEL YAMIN 
(Manuscript received March 7, 1972) 


A computer program is described which generates every intersection 
jigure resulting from the superposition of two closed polygon-like plane 
figures, each consisting of an arbitrary number of line segments or circular 
arc segments. Each intersection figure 1s assigned to one of four regions 
of the plane, representing the union, the intersection, and the two “exclusive- 
oR’s” formed by the pair of input figures. The two input figures may 
intersect or be tangent at any number of points and may have sections 
of coincident boundaries. No grid approximation is used. The program 
operates in two stages: the first stage analytically finds and classifies every 
point of intersection or tangency of the figures; the second stage regards 
these points as the nodes of a graph and applies an algorithm which causes 
each intersection figure to be traced just once. 


I. INTRODUCTION 


In the course of a project related to computer-aided integrated 
circuit mask design, it became necessary to describe the configuration 
which is formed when two polygon-like plane figures are superimposed 
on one another. Two closed figures in a plane divide the plane into four 
regions: inside both figures, inside the first but not the second, inside 
the second but not the first, and inside neither. Each region may consist 
of one or more figures. Assuming that the original figures consist of 
an arbitrary number of sides, each of which may be a line or circular 
arc segment, it was desired to describe every figure resulting from their 
intersection and to assign each to one of the four regions. 

The simple approach of establishing a grid of points and determining 
which sets of points are included in each region was not considered 
applicable, because the number of grid points required for sufficient 
resolution would have led to excessive computation. Instead, the 
following approach was used. Pairs of sides, one from each input figure, 
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were investigated analytically for points of intersection or tangency, 
or colinear sections. From the resulting information, a table of inter- 
section points was developed. Included in this table was information 
concerning the type of the intersection: crossing, tangency, beginning 
or end of a colinear section resulting in a crossing, or beginning or end 
of a colinear tangency. With the intersection table available, the col- 
lection of line and are segments constituting the superimposed figures 
could be considered as a graph, the intersections being the nodes. 
An algorithm for tracing the graph so as to generate every interséction 
figure once and only once was developed. The tracing paths thus 
determined were used to select appropriate analytical data from the 
original figure descriptions so as to generate descriptions of the new 
figures in the same numerical formats as the original figure descriptions. 
Classification of the figures was an automatic result of the tracing 
algorithm. 


II. NUMERICAL REPRESENTATION OF FIGURES 


Both input and output figures are described as sequences of sides 
in a numerical format designed for efficiency of computation rather 
than compactness of storage. The rotational tracing sense (clockwise 
or counterclockwise) of each figure is specified, and thus each side has 
a direction. A side, which may be a line or circular are segment, is 
represented by eight sequential words in memory: the starting co- 
ordinate pair, parameters of the analytical equation, and the terminal 
coordinate pair. Since the terminal point of one side is the starting 
point of the next, each additional side requires six additional words 
of storage. Each figure is explicitly closed with a side, which terminates 
at the starting point of the first side. 

The side sequence of each figure is preceded by a header containing 
information about the figure as a whole: its extreme z and y coordinates, 
allowing for projecting arcs; its area; the number of sides; and the 
sense (cw or ccw) of the entire figure. The header also contains pointers 
which can be used to associate the figure with externally tabulated 
information or can be used to associate it with other figures in one 
or more linked lists. Figure descriptions are stacked sequentially in 
memory; specific figures are located by pointers into this data structure. 
The exact arrangement of pointers and external tables depends on 
the application. 

Two limitations are imposed on the input figures; they must nowhere 
intersect themselves; and no side should be colinear with an adjacent 
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side. It is possible to override the second limitation with a preprocessing 
program which combines such sides and condenses the figure description 
appropriately. 


III. DEVELOPMENT OF INTERSECTION TABLE 


Given two figures described in the above format, it is desired to list 
all their points of intersection. This list of intersections will include 
the coorditiates of each intersection point, indices assigning it to a 
specific side of each of the two figures, and a code defining it as a cross- 
ing point, a tangency, one end of a colinear section topologically equiv- 
alent to a crossing, or one end of a colinear section topologically equiv- 
alent to a tangency. . 

It is typical of this entire analysis that much of the work is involved 
with the treatment of special cases. When the middle of one side is 
intersected by the middle of another, it is no trouble to make the 
appropriate entry in the intersection list. The problems arise when 
intersections involve more than one side of each of the input figures, 
as when they occur at corners, or when colinear sections ‘‘wrap around” 
several sides of each figure. To handle these cases, the intersection 
table must be developed in stages, entries ambiguous at a given stage 
being flagged as such and clarified in the next stage. 

One of the input figures is arbitrarily declared the operating, or A, 
figure; the other the passive, or B figure. The elementary process of 
the intersection table development is the solution of the side equation 
of one side of A with one of B to find any intersections, tangencies, or 
colinearities within the bounds of both segments. A solution which is 
detected just at the end of one of the sides can not generally be classified 
as a crossing or tangency at this stage. The second stage is to repeat 
the process, with the same side of A, for every side of B, so that a table 
of points of intersection for this one line or arc segment with B is created. 
At this stage, while ambiguities resulting from intersection at a corner 
of B have generally been cleared up, solutions at the ends of the A test 
side are still not finally classified. The third and last stage is to repeat 
the above for every side of A, pairing up unclassified solutions which 
correspond to intersections at corners of the A figure. 

Each of these stages of analysis is performed by a subroutine which 
loops on the previous one. Thus, a subroutine called INSECT returns 
the point or points of intersection of two line or are segments (that is, 
those solutions of their equations which lie within the bounds of both 
segments) or reports their colinear coincidence. SLASH returns the 
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points of intersection of a line or arc segment with a closed figure 
consisting of such segments; classifying as crossings, tangencies, etc., 
those intersections which do not involve the ends of the segment and 
returning special codes to provide what information it can about the 
others. INTLST returns the points of intersection of two closed figures 
consisting of line and are segments, each intersection being classified 
as a crossing, a tangency, or the start or end of a colinear tangency 
or crossing. This is the “intersection list’? previously described; the 
intersections appear in the order in which they would be encountered 
when tracing the A figure. Two special situations are recognized: if 
the figures nowhere intersect, an intersection list of zero length is 
returned; if the figures are coincident (everywhere colinear), a special 
code is returned. 

In principle, this procedure requires a number of calls to INSECT 
equal to the product of the number of sides of the two figures. However, 
each of the subroutines applies a quick preliminary test for overlap 
of the circumscribed rectangles of the input line segments or figures. 
In most cases, this avoids a great deal of detailed computation. 

INTLST and SLASH use a complicated sieve of analytical tests 
to determine the nature of intersections which occur at figure corners 
and involve three, sometimes four, are or line segments. Typically, 
these tests must decide whether the intersections represent a crossing 
or tangency of the two figures, and the geometrical question to be 
resolved is usually whether the direction of incidence of one line or 
are segment on another segment, or on the point of intersection of two 
other segments, is from the left or from the right. In the case of lines, 
a vector cross product is used to make this determination; in the case 
of ares, which may be tangent at a point of intersection, logic involving 
the direction of the arc center is required. 

All three levels of geometric analysis make use of a quantity called 
the “error margin,” which is necessary to allow for computational 
inaccuracies. The assumption is made that two points closer together 
than this small distance are in reality the same point. The error margin 
differs from a grid approximation in that its magnitude has no influence 
on computation time. 

One more item of information must be added to the intersection 
list for the purposes of figure tracing. This is the incidence direction 
of the A figure on the B figure, and of B on A, at every intersection. 
The incidence is a two-valued parameter; one value corresponds to 
entry of one figure into the other at a crossing intersection or inside 
incidence at a tangency; the other value corresponds to exit at a crossing 
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or outside incidence at a tangency. There are obviously four types of 
intersection, corresponding to the combinations of the two possible 
incidence values of A-on-B with those of B-on-A. It is convenient to 
determine these values as part of the program next to be described, 
rather than as part of INTLST. 


IV. GENERATION OF INTERSECTION FIGURES 


The development of the intersection table described above has been 
a geometrical problem in that finding the coordinates and the nature 
of each intersection has required knowledge of analytical parameters 
describing the exact position of each geometrical entity in the plane. 
With the intersection table in hand, the problem of tracing each inter- 
section figure becomes topological in nature. The overlapping polygons 
may be considered as a graph, that is, a set of intersection points, 
or nodes, between some of which exist connecting paths, or edges. 
For the purpose of figure-tracing, it is immaterial whether an edge 
between two intersections consists of one or many sides, or whether 
they are lines or arcs. It is enough to know to which figure the edge 
belongs. 

An algorithm has been devised by which each circuit of the graph 
which corresponds to an intersection figure may be traced once and 
only once. The sequence of edges constituting such a circuit may then 
be used to select analytical information from the input figure descrip- 
tions, from which a description in similar format of the intersection 
figure represented by the circuit may be assembled. Since each edge 
contains no intersections, it must belong to a single figure, which is 
known. Its start and end points, and the figure sides upon which they 
lie, are associated with its terminating nodes as listed in the intersection 
table. The tracing direction of a given edge may be the same as or 
opposite to the direction of the original figure. The generation of the 
intersection figures then consists essentially of copying sequences of 
sides or partial sides, forward or backward, as guided by the edge 
data and the edge sequence in the specified circuit. The region of the 
plane to which the intersection figure belongs is determined from the 
incidence information. 

The most direct way to implement such a procedure is based on the 
fact that each edge is part of two intersection figures. The method is 
to list, in cyclic order, all the edges radiating from each node. Then, 
starting from any node along an edge, proceed to the next node and 
transfer to the adjacent edge in cyclic order, clockwise or counter- 
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clockwise as appropriate, so that the trace remains within a “‘cell” 
of the graph. Continue until the starting node is reached, then repeat 
with a new starting node. The complete set of intersection figures will 
be traced when each edge has been traversed just twice (a colinear 
section being considered one edge). Unfortunately, no explicit numerical 
description of the edges connecting the nodes is available at this stage 
in the analysis; the necessary information must be computed from the 
input figure descriptions, the side blocks of which do not, in general, 
explicitly represent edges of the graph. 

The algorithm actually used, therefore, centers its attention on the 
nodes (intersection table entries), rather than the edges. As the trace 
encounters each node, the algorithm specifies the next node which 
must be encountered to keep the trace within a “cell” of the graph; 
the numerical description of the connecting edge is only then generated 
and posted directly to the output figure description. In principle, each 
node must be encountered four times for completeness, once for each 
of the figures of which it is a corner. The implementation of this test 
for completeness is complicated by special considerations which arise 
at nodes which are tangencies or colinearities; while the verbal descrip- 
tion of these special cases in Section 4.1 appears complex, their im- 
plementation in computer code requires only minor modifications in 
the program flow. 

The circuit tracing and intersection figure generation operations 
are embodied in a single subroutine called OVRLAP. This subroutine 
calls INTLST to establish the intersection table, in which the inter- 
sections are listed in A-figure order. For efficiency, it is desirable to 
have a table of intersections in B-figure order. This could be created 
by another call to INTLST, but it is much more efficient to produce 
it by reordering a copy of the A-figure intersection table. The two 
tables are crossreferenced so that any entry in one can be located in 
the other. The OVRLAP subroutine also finds, and tabulates, the 
A-on-B and B-on-A incidence values at each intersection as described 
previously. 

OVRLAP is a modular program, with modules (procedures) which 
perform elementary functions such as: 


(2) Extract a side from the A or B figure. 
(22) Insert a side in an output figure description. 
(112) Replace a side description by a description of the same side, 
but backwards. 
(tv) Given a side description and a point on the side, truncate the 
side to return the first or second segment. 
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Other modules use these to perform higher-level functions, such as: 

(zt) Insert in an output figure description a sequence of sides from 

one intersection to another, derived from the A or B figure as 
required, forward or backward as required. 

(772) Generate an output figure by tracing a complete circuit, ac- 
cording to the circuit-tracing algorithm, from a given inter- 
section back to itself. 

(417) Copy an input figure bodily as an output figure, making ap- 
propriate changes in the figure header. 


4.1 Circuit Tracing Algorithm 


First, for simplicity, consider a pair of intersecting figures in which 
all the intersections are point crossings. Tangencies and colinearities 
introduce complications which require modified treatment. 

(2) Starting at any intersection, trace two circuits. One is initiated 
by tracing forward on the B figure, the other by tracing backward 
on the B figure. Before starting these traces, compute a circuit incidence 
code for each. This code has four possible values, one for each of the 
possible combinations of A-on-B and B-on-A incidence values. For 
the B-forward trace, it combines the two incidence values (B-on-A 
and A-on-B) stored for the starting intersection; for the B-backward 
trace, it combines the A-on-B value with the opposite of the B-on-A 
value. 

The procedure of “‘tracing’’ consists of going to the next intersection 
on the appropriate figure in the chosen direction. The ‘‘next inter- 
section” is an adjacent item in the A-figure or B-figure intersection 
table, as appropriate. The last item on each list is considered ‘‘adjacent’’ 
to the first. 


(21) At the next intersection: 


If tracing on B, switch over to A backwards. 

If tracing on A, switch back to B in the starting direction. 

Label every intersection encountered in this way with the circuit 
incidence code for the current trace. Each item block in the B-figure 
intersection table has space to store four such code labels. 

Repeat step (zz) until the starting intersection is encountered, com- 
pleting the circuit. 

(222) On completing both circuits from the current starting inter- 
section, go to the next intersection on the A figure and proceed as 
under (2) above. However, do not perform any trace for which the 
starting intersection is found already to be labelled with its circuit 
incidence code. 
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(iv) Continue this procedure, making each intersection along the 
A figure the starting intersection in turn, until every intersection has 
been the starting intersection. 


On termination of this procedure, every circuit on the graph which 
corresponds to an intersection figure will have been traced just once. 
Moreover, the circuit incidence code associated with each figure directly 
assigns the figure to one of the four regions into which the two inter- 
secting figures divide the plane, thus: 


Circuit Incidence Code Means Region of Plane 
A enters B, B enters A Inside A, inside B (A and B) 
A enters B, B leaves A Outside A, inside B (B and not-A) 
A leaves B, B enters A Inside A, outside B (A and not-B) 
A leaves B, B leaves A Outside A, outside B (not-A and 


not-B) (A or B) 


(The last category represents the logical union of the two figures. 
Such a figure may be a “hole” in the union, or the outer boundary of 
both figures. In the latter case, the inside of the figure is logically the 
complement of not-A and not-B, namely, A or B). 

This algorithm follows from the fact that each crossing intersection 
forms a corner of four figures, each in a different region of the plane 
as described above. Step (7) ‘‘stakes a claim” on two of these figures. 
(The other two will be or have been claimed by the previous A inter- 
section.) Step (zz) assures that the trace stays within the “claimed” 
figure. Step (277) moves the procedure forward, assuring that no figure 
previously ‘“‘claimed”’ is retraced. Since every intersection figure is 
“claimed” by at least one intersection, step (7v) assures that every 
figure is eventually described. 

A point of tangency between two figures differs from a crossing inter- 
section in that two of the four figures which impinge on it belong to 
the same region of the plane. The following special considerations 
apply to ta igencies: 


(z) On encountering a tangency when tracing on the A figure, 
switch over to the B figure only if the B-on-A incidence value at the 
tangency is the same as at the starting intersection for a B-forward 
trace or opposite for a B-reverse trace. When tracing on the B figure, 
switch to the A figure only if the A-on-B incidence is the same as at 
the starting intersection. If this condition is not met, go on to the 
next intersection. 
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(it) When switchover from B to A is made at a tangency, do not 
label the tangency with the circuit incidence code. 

(cit) When it is the turn of a tangency to be the starting intersection, 
only one tracing circuit, at most, rather than two, is originated. The 
permissible tracing direction of B is determined by the following ‘‘inside- 
outside” rule: If the rotation senses of the two input figures are the 
same (both ew or both cew), trace backward on B if the B on A and 
A on B incidence values at the tangency are the same, and forward 
on B if they are opposite. If the rotation senses are opposite, the rule 
is reversed. Of course, just as with crossing intersections, no trace is 
performed if the tangency is already labelled with the appropriate 
starting incidence code. 


Hither a tangency or a crossing of the two input figures may occur 
as a colinear section. Such a colinear section appears in the intersection 
tables as two adjacent intersections, one labelled as the “‘start’’ and 
the other as the “end” of the colinearity. ‘Start’ and “end” refer to 
the direction sense of the A figure. The following special considerations 
apply to colinearities. 


(t) On encountering either end of a colinearity while tracing a 
circuit, take the appropriate switchover action immediately, following 
the rules described above for point crossings or point tangencies as 
appropriate. 

(21) Never begin tracing a circuit at the “start” of a colinearity; 
go to the ‘‘end”’ of the colinear section and make that the starting 
intersection as described above. 


(a) (b) 


Fig. 1—Two input figures: (a) composed of 16 sides, 12 of them lines and 4 arc 
segments; (b) composed of 15 sides, 11 of them lines and 4 arcs. ; 


1604 THE BELL SYSTEM TECHNICAL JOURNAL, SEPTEMBER 1972 


Fig. 2—Graph resulting from superposition of Figs. la and 1b. 


(it) If a colinear region is being traced, having been entered at 
the “end” point (whether by switchover or by initiation of circuit 
tracing), the next intersection encountered will be the “‘start’’ point. 
Do not switch figures at this point; just go on to the next intersection. 
Any switchover required will already have been done at the “end” 
intersection. 


4.2 Generation of Output Figure Description 


The algorithm just described sets out a procedure for stepping from 
node to node in the intersection lists in such a way as to describe certain 
unique circuits. In essence, it tells one to choose an appropriate item 
in the A-figure list; find the corresponding entry in the B-figure list; 
follow the B list item by item, in a specified direction, until certain 
conditions are met; switch to the corresponding entry in the A-figure 
list; follow this list backward, item by item, until certain conditions 
are met; switch to the corresponding entry in the B-figure list; and 
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continue switching back and forth in this way until the starting point 
is attained. 

The numerical description of the output figure which is the objective 
of this program is generated simultaneously with the execution of this 
algorithm. As was mentioned earlier, the program contains a module 
which, given two specified intersection nodes, will generate a sequence 
of side descriptions representing a path from the first intersection to 
the second, on the A or B figure as specified, and forward or backward 
as specified. This module is invoked every time a switch from one 
intersection list to the other is made. The two intersections specified 
are that at which the list was entered, and that from which the switch 
is being made (if there are no tangencies or colinearities, these are 
adjacent intersections). The figure is that represented by the inter- 
section list, and the direction is that in which the list was followed. 
The resulting side sequences are stacked one after another in memory, 


Fig. 3—Figures classified as ‘‘A or B’’-the union of two input figures. 
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and when the circuit trace is complete, so is the intersection figure 
description. . . 

An anomaly can occur as a result of the circuit-tracing algorithm: 
two successive sides of an intersection figure may be returned which 
are colinear extensions of one another. It is easy to recognize this 
situation, and the program module which transfers individual sides 
to the output figure simply combines any two such sides. 

Before a circuit trace is initiated, the intersection figure description 
must be initialized by allocation of space for its header. After the circuit 
trace is complete, certain header items must be filled in, particularly 
extreme coordinate values and the area. The area of a figure is found 
by a sum-of-trapezoids technique, with corrections for ares. The area 
returns as a signed number, the sign denoting the rotational sense 
(ew or cew) of the figure. The header includes a code which specifies 
to which of the four regions of the plane the figure belongs. 

It can occur that two input figures have no crossings at all, but 
only tangencies. One may be inside the other, or they may be external 
to each other. In this case, the input figures are themselves intersection 
figures, but they will not be traced by the above algorithm. Therefore, 
they must be copied bodily as output figures and assigned by simple 
logic to an appropriate region of the plane. The same is true for figures 
which do not touch each other at all. 

Another problem arises when two figures have one point tangency 
as their sole intersection. A figure will be generated which is wasp- 
waisted or self-tangent, either of which violates the rule against self- 
intersection. The program returns these descriptions regardless, and 
the calling program which uses OVRLAP must watch out for this 
situation. Finally, when one figure lies within another without any 
contact, one region of the plane (inside one figure but outside the 
other) cannot be described as a simply connected figure, though its 
area is easily calculated. 


4.3 Organization of Output Data Structure 


Intersection figures generated as described in the previous section 
are stored sequentially in a region of memory allocated for output. 
To be useful, this data structure must be organized and indexed in 
some way. The following system has been used: 


(t) All the figures of a given type (that is, belonging to the same 
one of the four regions of the plane) are chained together in a 
linked list by a pointer in the figure headers. 
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(2t) A summary table is provided. This table has four blocks, each 
corresponding to one of the four regions of the plane. Each 
block contains: the number of figures of this type, the total 
area of these figures, and a pointer to the first figure in the 
corresponding linked list. 

(i772) Each figure header has a pointer to the appropriate block of 
the summary table, and one to the next sequential figure descrip- 
tion in the data structure. 


V. DETAILS OF IMPLEMENTATION 


The program described in this paper has been coded in FORTRAN IV 
and compiled and executed on the General Electric 635 computer. 
It has been coded as a subroutine, called OVRLAP, with three ar- 
guments: The locations of the two input figure descriptions and the 
starting location of a region of storage available for output. In addition, 


Fig. 4—Two regions of overlap, ‘‘A and B.” 
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common blocks must be provided for the summary table described 
above, for a summary table (of different form) which controls the data 
structure in which the input figures are located, and for the error 
margin quantity discussed in Section III. The common block con- 
taining the summary table has three additional locations. One reports 
back the total number of figures generated; another carries a code for 
the mode of intersection (no contact, tangencies only, coincident figures, 
or normal intersection). The third is set on input as the total number 
of locations available for output figure descriptions and returns as a 
pointer to the first unfilled location. An error return occurs if the 
output data structure tries to overflow the space allocated for it. 
Figures 1 through 6 are examples of the operation of this program. 
Figure 1 shows the two input figures. The first has 16 sides, 12 of them 
lines and the other 4 arc segments; the second has 15 sides, 11 of them 
lines and the other 4 ares. Figure 2 shows the graph resulting from 
their superposition. There are 15 points of intersection, and every 


Fig. 5—Four regions of A not covered by B—A and not-B. 
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Fig. 6—Five regions of B not covered by A—B and not-A. 


type is represented: point crossings, point tangencies, colinear crossings, 
and a colinear tangency. Thirteen intersection figures, including the 
outer boundary of the pair, are generated. Figure 3 shows the figures 
classified as ‘“‘A or B” (the union of the two input figures). Note that 
there is an inner figure, or hole, which can also be described as bounding 
a region “not-A and not-B”. Figure 4 shows two regions of overlap, 
“A and B” (the intersection, in the logical sense, of the two figures). 
Figure 5 shows four regions of A not covered by B (A and not-B) and 
Fig. 6 shows five regions of B not covered by A (B and not-A). 

Execution of this analysis by the OVRLAP subroutine required 
0.75 second of processor time on the GE 635, the memory cycle of which 
is one microsecond. This does not include time spent on physical input, 
generation of the input data structure, or physical output. 

OVRLAP and all its required subroutines and common blocks, 
exclusive of system library subroutines, occupies about 11000 (decimal) 
words of storage. 


VI. SUMMARY 


A computer program has been developed which, given two polygon- 
like figures each consisting of an arbitrary number of line segments 
or circular arc segments, generates every intersection figure resulting 
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from their superposition and identifies its logical relationship with the 
input figures. The two input figures may intersect or be tangent at 
any number of points and may have sections of coincident boundaries. 
Each input figure is described as a sequence of bounded analytical 
expressions and the output figures are generated in the same numerical 
format. No grid approximation is used; the program operates by ana- 
lytically finding every point of intersection or tangency of the figures; 
these points are the nodes of a graph through which the program traces 
according to an algorithm which causes every intersection figure be 
traced just once. A FORTRAN IV implementation of this program 
took 0.75 second on a GE 635 computer to resolve the superposition 
of two figures, of 15 and 16 sides respectively, into 13 resultant figures. 
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A computer simulation of two-dimensional mobile radio systems arranged 
with square coverage areas on a square grid and using dynamic channel 
assignment techniques is described. Parameters for the simulation are: 
(t) 729 distinct coverage areas (27 on a side), (i7) 160 radio channels 
and (at) a radio channel reuse interval of every fourth coverage area. 
Three different channel assignment strategies are considered and the 
results are compared to previous one-dimensional simulated systems and 
to a fixed channel assignment system. At call blocking rates below about 
10 percent, the two-dimensional dynamic channel assignment systems carry 
more traffic and produce fewer forced call terminations at coverage cell 
boundaries than do fixed channel assignment systems. For example, at a 
blocking rate of 1 percent, the traffic carried, TC, expressed in Evrlangs 
per channel per coverage area by the various systems are as follows: fixed 
channel assignment systems, TC = 0.44, one-dimensional dynamic channel 
assignment system, TC = 0.62, two-dimensional dynamic channel assign- 
ment system, TC = 0.68. 


I. INTRODUCTION 


When a large pool of radio channels is available for use in a large-scale 
multiple base station mobile radio system, its assignment to individual 
base stations on the basis of instantaneous demand can improve systems 
performance.” ° Previous studies of coordinated dynamic channel 
assignment in one-dimensional systems (that may be used along an 
expressway or major air route for example) have shown that at blocking 
rates below 10 percent these assignment procedures produce increased 
channel occupancy over systems using the same channels but with 
channel allocation to each base station fixed. 

The most general form of dynamic channel assignment assumes that 
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any channel can be used at any base station. Efficient use of channels 
requires the simultaneous assignment (reuse) of the radio channel in 
radio coverage areas which are spaced as close together as possible 
without incurring excessive cochannel interference. Mobile radio 
systems which reuse channels within a metropolitan area often are 
referred to as small-cell or small coverage area systems.* ° 

The analysis of dynamic channel assignment systems which employ 
complex channel assignment algorithms appears to be intractable. lor 
this reason, the performance of these systems is most readily determined 
by large-scale computer simulation. The overall system performance 
of two-dimensional dynamic channel assignment mobile radio systems 
is expected to be similar to that of the one-dimensional systems pre- 
viously studied.’'? However, the performance of the one-dimensional 
systems was shown to depend upon the total number of channels 
available to the radio system as well as the average number of channels 
available in each coverage area. Thus, since the channel reuse situation 
is more complex for two-dimensional systems, some difference in the 
performance of the two-dimensional systems was expected. 

This paper describes a computer simulation of a large-scale two- 
dimensional mobile radio system using different dynamic channel 
assignment strategies and compares the performance characteristics 
of these systems with both a one-dimensional dynamic channel assign- 
ment system and with a fixed channel assignment system. 


II. THE SIMULATION MODEL 


2.1 The Model of the System Configuration 


The system simulated consists of a set of square radio coverage areas 
arranged to completely cover a square grid with no overlap as illustrated 
in Fig. 1. For this model, a channel used in a specified coverage area, 
such as the shaded area at coordinates m, n may be reused anywhere 
on or outside of a specified square ring which surrounds that particular 
coverage area. This reuse ring is defined by those coverage areas which 
have either an x or y coordinate which is a specified integral number of 
coverage areas separated from the center coverage area. This specified 
interval of coverage areas is referred to as the channel reuse interval. 
In the example in Fig. 1, the cross-hatched coverage areas on the reuse 
ring surrounding the coverage area at x = m, y = n are at a reuse 
interval of 4. Several coverage areas on the reuse ring may also use the 
same channel as long as all coverage areas are separated by at least 
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Fig. 1—Coverage area configuration for the simulation. 


the reuse interval. Thus, in the example, the maximum number which 
can use the same channel simultaneously on the reuse ring is eight. 
lor the simulation the coverage squares were specified to be 2 miles 
on a side. 

The “local” mean (average over several wavelengths) signal strength 
in an urban area measured on a circle of constant radius from an omni- 
directional base station has a log normal distribution.’** The variance 
of the log normal distribution is sufficiently large (¢ = 5 to 9 dB) so 
that circular coverage areas are seldom realized. Since traffic handling 
capabilities of dynamic channel assignment systems are the primary 
concern of this study, it is desirable to define geometric coverage areas 
which completely cover a region with no overlap. The only candidates 
are triangles, rectangles, and hexagons. Of these, the geometry asso- 
ciated with squares is simplest to manipulate. Since the propagation is 
highly irregular, radio channels obviously may not be used simul- 
taneously in immediately adjacent coverage areas regardless of the 
geometry defining the coverage area. Thus, square coverage areas on 
a square grid with the reuse intervals described previously provide a 
reasonable model for these traffic studies. Such a system could be 
approximated by placing base stations at the corners of the coverage 
areas and radiating into them with directional antennas. 
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2.2 The Stimulus to the Simulated System 


The initiation of call attempts and the movement of vehicles in this 
simulation is almost identical with that described previously’ for a 
one-dimensional simulation. Call attempts are generated in the simu- 
lation as a Poisson process in time which has a controllable average 
attempt rate in each coverage area. For this study, initial call attempts 
from vehicles are uniformly distributed in both the x and y coordinates 
and the x and y locations are independent. Call attempts are tapered 
at the edges of the simulated universe, and the performance data 
accumulated only from the central portion of the system, so that the 
results are representative of an infinite system. 

In this method of simulation, vehicles making calls are identifiable 
entities the locations and movements of which are stored in the com- 
puter. Vehicle velocities (speeds and direction) are chosen randomly 
from a population having statistical characteristics which can be 
prescribed. Vehicle motions are primarily parallel to one of the coordi- 
nate axes with 45 percent of the vehicles having an x velocity component 
only, and 45 percent a y component only. Ten percent of the vehicles 
have both x and y velocity components which are mutually independent. 
The velocity components (x and y) have a truncated Gaussian distri- 
bution with a zero mean, a standard deviation of 30 miles per hour, and 
a maximum velocity magnitude of 60 miles/hour (see Figure 8 of Ref. 1). 

Call attempts that are assigned a channel remain on in the system 
for call durations which are taken from a random population with a 
specified distribution. Part of the data to be presented later is for a call 
duration distribution which is exponential with a mean of about 98 
seconds. Other data are for a call duration distribution which is a 
truncated Gaussian with a minimum call duration of 30 seconds, a 
maximum call duration of 10 minutes, and a true mean of 103.5 seconds. 
The mode of this distribution is 90 seconds (see ligure 7 of Ref. 1). 
Some of the vehicles which cross coverage area boundaries have their 
calls prematurely terminated by the system because of the unavailability 
of channels in the new coverage area. This effect shortens some calls and 
thus perturbs the specified distribution slightly. This is discussed in 
more detail later. 


2.3 The Simulated Operating Systems 


The method of handling call attempts, calls in progress and call 
terminations is similar to that described in Ref. 1 for a one-dimensional 
system. The first step in processing a new call attempt is to assign a 
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coverage area. Since coverage areas are geometrically defined and the 
vehicle coordinates known, the appropriate coverage area is readily 
determined. In the simulation, vehicles within a particular coverage 
area are required to be served by channels available in that coverage 
area. The next step is a channel search which tries to find a radio channel 
to serve the awaiting call attempt. Some channel search procedures 
are discussed in detail later. If no channel can be found, the call is 
immediately blocked and cleared from the system. This blocking 
strategy is the same as that used in deriving the Erlang B (block call 
cleared) telephone traffic formula.” 

The simulated system detects vehicles crossing coverage area bound- 
aries and checks to see whether the call can continue on the original 
channel or on a new channel. If no channel is available, the call is 
immediately forced to terminate and is cleared from the system. In an 
actual system, the crossing of coverage area boundaries will not be 
critical and, hence, some delay could be tolerated in finding a new 
channel for boundary crossing vehicles. When a significant number of 
the calls handled by the system experience a boundary crossing, the calls 
forced to terminate have an effect on other system performance param- 
eters. With the coverage area size and the velocities used in this simula- 
tion the effect of boundary crossings on traffic carried and call blocking 
is small. 


2.4 Channel Assigninent Strategies 


Information which identifies the coverage area and the channel being 
used by each active subscriber in the network is stored in the memory 
of the simulated system control computer. In addition, a list is kept of 
every channel being used in each coverage area. The second list, although 
redundant, is in a form which permits rapid access to the information 
needed for the dynamic channel assignment procedures. A channel 
search is initiated after the coverage area of a call attempt has been 
determined. 

The procedures used in the channel assignment strategies are best 
described in terms of an example which can be referred to Fig. 1. Assume 
that a call attempt located at coverage area v = m; y = n is awaiting 
channel assignment. The first step is to determine which channels 
Gif any) are available to serve the call attempt. This is accomplished 
by compiling a list of the channels which are not being used currently 
within any coverage area surrounding the designated coverage area 
out to but not including the reuse ring defined previously. If more than 
one radio channel is available within this reuse interval, then some 
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channel assignment strategy must be applied to determine the channel 
that should be assigned. The main object of a channel assignment 
strategy is to increase the channel occupancy or to optimize some 
other system parameter. Briefly, the two-dimensional channel assign- 
ment strategies considered here are as follows: 


First Available 


This strategy assigns the first available radio channel encountered 
during a channel search. This strategy obviously minimizes the system 
computation time and as shown later its performance may be adequate 
for some applications. 


Selection of a Channel With Maximum Usage on the Reuse Ring (RING) 


A search is made through the list of available channels to determine 
which of these channels is currently in use in the most coverage areas 
lying on the reuse ring. If more than one channel has this maximum 
usage, an arbitrary selection of one of the channels is made to serve the 
call attempt. If none of the available channels is in use on the reuse 
ring (an infrequent event as illustrated in Tig. 8), then selection is 
made on a first available basis. This strategy is described in detail 
in Ref. 10. 


Selection of a Channel With Maximum Usage on a Side(s) of the Reuse 
Ring (Orthogonal NN) 


A search is made through the list of available channels to determine 
which channel(s) meets the following criteria: 


(2) On some side of the square reuse ring channel usage is a maximum. 

(iz) Of those channels having a maximum usage on some side, channel 
usage 1S again a maximum on a side adjacent (orthogonal) to that 
first maximum usage side. 


If none of the available channels are in use on the reuse ring, channel 
assignment is made on a first available basis. If only one channel has 
a maximum channel usage on some side, that channel is assigned 
without regard to the usage on other sides. If more than one channel 
has maximum usage on some side and equal secondary maximum usage 
on an adjacent side, one of these channels is selected at random. For 
this simulation optimization was not carried to the third and fourth 
sides. This strategy is one possible extension of the best performing 
one-dimensional channel assignment strategy” previously described. 

Figure 2 illustrates some of the differences in these channel assignment 
strategies. In this figure, assume that channels are designated by the 
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letters A through E. Channel A is not available for assignment because 
it is being used in coverage areas within the reuse interval. For channel 
selection by the RING strategy, channel B is used four times on the 
reuse ring, C six times, D five times, and E four times. Thus the RING 
strategy selects channel C for assignment. 

For the orthogonal NN strategy the maximum side usage is three, 
and occurs for channel D on the upper side and channel E on the left 
side. Channel usage on a corner is included in the usage count for each 
side which includes that corner. for example, the channel D usage on 
the upper left-hand corner of Fig. 2 results in a channel usage count of 
three for the upper side and a count of two for the left side. Thus, the 
adjacent side usages for channel D are two on the left and one on the 
right. The corresponding adjacent side usages for channel E are only 
one each on the top and bottom. Therefore, channel D has the maximum 
adjacent side channel usage and is thus selected for assignment by the 
orthogonal NN strategy. 

If the channel assignment search is assumed to proceed in alphabetical 
order, the first available strategy would select channel B for assignment. 


2.5 System Parameters 


The simulated operating systems were run for a 27 by 27 grid of 
coverage areas which were 2 miles on a side. The systems had available 
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160 duplex radio channels. With the assumed reuse interval of 4, this 
results in 10 channels (160/16) available on the average at each base 
station. Of course, depending upon the system activity and the instan- 
taneous demand, some coverage areas will have more than 10 channels 
in use simultaneously and some will have fewer. Because of the uniform 
spatial and temporal statistics of the call-attempt process, averages of 
system parameters taken over a large number of coverage areas at a 
particular time must yield the same results as averages taken over an 
extended time interval for one coverage area. Thus, the stability of the 
statistics for performance parameters depends upon the product of 
the number of coverage areas included and the run time. Statistics 
were taken from only the central 225 coverage areas (15 by 15) to avoid 
any effects from the edges of the finite system. Call loading and system 
response at the edges of the finite system were tapered in a way which 
insured that these central areas were operating as members of an 
infinite set. The simulation was started initially with no calls on in the 
system and statistics of parameters were monitored until the initial 
transients died out. The time required for stability was about 2 minutes. 
Data were collected after stabilization for about 7 minutes (in simulated 
time). This run time and spatial size was sufficient produced data 
points which lie on smooth curves and have very little scatter. 

The data presented in the following sections were obtained by 
counting and storing the number of actual events which occurred in 
the stabilized system. Some of the events counted were new call- 
attempts, new call-attempts blocked, coverage boundary crossings and 
calls forced to terminate after crossing coverage boundaries. In addition, 
the actual number of calls in progress in the system was counted and 
stored periodically during each simulation run. 


II. COMPARISONS OF SYSTEMS PERFORMANCE 


The computer simulation was run at several call attempt rates for 
the channel assignment strategies and system parameters described in 
Section II. The data from the simulation, which show directly the 
relationship between channel occupancy (traffic carried per channel 
per coverage area) and the blocking of new call attempts, are plotted 
in Fig. 8. The solid data points are for the exponential call duration 
distribution and the open data points are for the truncated Gaussian 
call duration distribution described in Section 2.2. The triangles and 
circles are for the two-dimensional systems with a reuse interval 
of four. The solid curve is drawn through the data points for both the 
RING and the orthogonal NN two-dimensional channel assignment 
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strategies. The performance of these two strategies in terms of these 
parameters is nearly identical, as indicated by the fact that the solid 
triangles at 1.1 percent and 6.6 percent blocking lie nearly on top of 
the solid circles. All of the points for these two strategies fall on the 
smooth curve over the blocking range indicated (1 through 30 percent). 
In fixed channel assignment systems for which the Erlang B telephone 
traffic formula is applicable, it has been shown" that the shape of the 
call duration (holding time) distribution has no effect on the relationship 
between traffic carried and the blocking of new call attempts. The 
open circle at 2 percent blocking that lies on the solid curve is for a 
different call duration distribution than the solid points as noted in 
the key. This indicates that for these dynamic channel assignment 
strategies the performance parameters are also not affected by the 
shape of the call duration distribution. The inverted triangles connected 
by the long dashed curve are for the two-dimensional first available 
channel assignment strategy. This strategy results in slightly lower 
channel occupancy for a given blocking. It is much simpler to implement, 
however, and may be acceptable in some applications. 

The short dashed curve through the open square data points is from 
a one-dimensional channel assignment system with a reuse interval 
of 4, as described in detail in Ref. 2 (NN strategy). The solid squares 
are for the same one-dimensional channel assignment strategy, but 
for an exponential call duration distribution. The fact that the points 
lie on the curve for the truncated Gaussian call duration distribution 
further illustrates the independence of the performance of dynamic 
channel assignment to the type of call duration distribution. Over a 
wide blocking range, the traffic carried by the two-dimensional system 
is more nearly constant than for the one-dimensional system. We would 
expect that at very low blocking rates the two-dimensional system 
would perform better (carry more traffic at a given blocking rate) 
because it has more channels available to the system (160 channels) 
than the one-dimensional system (40 channels) even though the average 
number of channels available per coverage area is the same. However, 
as seen from Iig. 3, the two-dimensional system does not perform as 
well as the one-dimensional system above a blocking rate of about 
1 percent. This is probably due to the fact that the simultaneous use 
of the same channel in coverage areas separated by exactly a reuse 
interval is more difficult to achieve in the two-dimensional area situation 
than on a one-dimensional line. To further check the performance of 
the simulations, the two-dimensional orthogonal NN strategy system 
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was run as a one-dimensional system by restricting the location of call 
attempts and vehicle motion to one axis. The resulting data points for 
this simulation were located along the performance curves from the 
previous one-dimensional simulation. 

The four data points denoted as ® are from a simulation of a one- 
dimensional fixed channel assignment system.’'” In this system, specific 
channels are allocated to each coverage area and the same channels are 
reallocated at coverage areas separated by exactly a reuse interval. 
For the simulated fixed channel assignment system, 10 channels were 
allocated to each coverage area and the reuse interval was 4. The 
performance of the fixed channel assignment system can be determined 
directly from the Erlang B telephone traffic equations.”’” The relation- 
ship between the traffic carried and call blocking from the Erlang B 
formula is the dotted line in Fig. 3. This curve represents the perform- 
ance of fixed channel assignment in both one and two dimensions since 
the curve depends only on the number of channels available per coverage 
~ area. The fixed channel assignment system performs worse at low block- 
ing than any of the dynamic channel assignment systems because it 
is less able to meet the peaking of call attempts in the randomly offered 
trafic. At high blocking rates, fixed channel assignment systems handle 
more traffic because channel reuse is fixed in an optimum configuration. 

The two-dimensional simulation using the RING strategy was also 
run for a reuse interval of five with 250 channels, and for a reuse interval 
of six with 360 channels. The resulting data points are indicated by the 
X and &) in Fig. 3. As was found in the one-dimensional case,” there 
was a small increase in the traffic carried at a specified blocking as the 
reuse interval increased even though the average number of channels 
available per coverage area was held constant. The factors affecting 
this increase in traffic carried are discussed in Ref. 2. 

Another measure of system performance is the behavior of the 
blocking and traffic-carried parameters plotted as functions of the new 
call-attempt rate as shown in Fig. 4. The symbols and the types of 
lines denote the same system parameters and channel assignment 
strategies as they did in Fig. 3. At low attempt rates, no new call- 
attempts are blocked in any of the systems and the traffic carried 
(average calls on per channel per coverage area) is a linear function 
of the new call-attempt rate. As the attempt rate increases, the fixed 
channel assignment system begins to block calls before any of the 
dynamic channel assignment systems. The two-dimensional systems 
are the last to experience significant (greater than 1 percent) blocking. 
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Fig. 4—Comparative performance of mobile radio systems. 


THE BELL SYSTEM TECHNICAL JOURNAL, SEPTEMBER 1972 


45 


40 


35 


30 


25 


20 


10 


1.2 


CALLS BLOCKED PERCENT OF NEW CALL ATTEMPTS 


DYNAMIC CHANNEL ASSIGNMENT 1623 


The curves labeled “calls on” illustrate that at high attempt rates the 
fixed channel assignment systems carry more traffic than the dynamic 
channel assignment systems. 

With no vehicle motion, the loading or traffic offered to these systems 
is the product of the new call-attempt rate and the mean call duration 
determined by the stimulus part of the simulation. One Erlang of 
offered traffic (one channel occupied 100 percent of the time) for the 
system unperturbed by motion is indicated on the call attempt axis. 
When vehicles cross coverage area boundaries and do not find channels 
available for continuing their calls, their premature termination reduces 
the actual average call duration. The fact that these systems must 
attempt to find new channels for some of these boundary crossing calls 
produces additional loading on the system. The actual average call 
duration of the exponential call duration distribution (98 seconds) is 
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Fig. 5—Forced call terminations relative to boundary crossings. 
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less than that for the truncated Gaussian distribution (103.5 seconds). 
This difference in call duration causes the actual loading on the systems 
to be different at a given new call-attempt rate and accounts for the 
displacement of the data points for the exponential call duration dis- 
tribution from those of the truncated Gaussian in Fig. 4. This effect 
is not seen in Jig. 3 because the traffic carried contains the actual call 
durations. The small difference in call durations resulted from different 
quantization effects in the methods used to generate samples from the 
two distributions. 

Since some of the boundary crossing calls in the dynamic channel 
assignment systems keep their originally assigned channels, the forced 
call termination rate is less than the blocking rate in these systems. In 
the fixed channel assignment system, all boundary crossing calls require 
a new channel, and thus the forced termination rate is the same as the 
blocking rate. In Tig. 5, the forced call termination rate expressed as 
the ratio of the number of calls forced to terminate in a time period to 
the total number of active call boundary crossings in that time period 
is plotted vs the blocking rate of new call-attempts. This function 
appears to be the same for all of the dynamic channel assignment systems 
included in Fig. 5. 

In Fig. 6, the ratio of the number of calls keeping the same channel 
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when crossing a coverage boundary to the total number of boundary 
crossings is plotted as a function of the new call-attempt rate. This 
relationship indicates the extent to which the originally assigned radio 
channel ‘‘floats’”’ with the vehicle. The behavior of the two-dimensional 
RING and orthogonal NN strategy are similar for these parameters 
also. Fewer callers keep their originally assigned channels when crossing 
boundaries in the two-dimensional systems than do in the one-dimen- 
sional systems. 


IV. OPERATING CHARACTERISTICS OF THE TWO-DIMENSIONAL DYNAMIC 
CHANNEL ASSIGNMENT SYSTEMS 


More than one channel must be available for assignment in a partic- 
ular coverage area if a channel assignment strategy is to improve 
system performance. Figure 7 shows the number of channels available 
for assignment expressed as the percentage of total channel assignments 
(new call attempts plus boundary crossing calls having channels re- 
assigned) attempted for several blocking rates. At low blocking rates, 
several channels are available to select from for most of the channel 
assignments attempted. As the loading on the system increases, fewer 
channels are available from which to choose. Thus, the effectiveness 
of the channel assignment strategies decreases as the system loading 
increases. Although slight differences exist in the histograms for the 
different channel assignment strategies, the overall availability of 
channels for each assignment strategy is about the same. 


4.1 RING Strategy 


The two-dimensional RING strategy attempts to maximize channel 
usage on the reuse ring. Figure 8 illustrates the degree to which this is 
achieved for two blocking rates. The percentage of total channel assign- 
ments made (new call attempts plus boundary crossings having channels 
reassigned) is given for each possible number of simultaneous channel 
usages on the reuse ring. The histograms of lig. 8 show that this dynamic 
channel assignment strategy is not able to maximize channel reuse as 
it serves the randomly-offered call attempts. It should be possible to 
reassign channels to calls in progress so that the channel reuse is more 
concentrated and thus the average channel occupancy is increased for 
a given blocking rate. 

The number of occurrences of channel assignments with each possible 
number of simultaneous channel usages on the reuse ring is expressed 


DYNAMIC CHANNEL ASSIGNMENT 1627 
30 


20 6.65% BLOCKING 





1.14% BLOCKING 


OCCURENCES PERCENT OF TOTAL 
ASSIGNMENTS MADE 





012 3 4 5 6 7 8 


NUMBER OF USAGE OF A CHANNEL ON REUSE 
RING FOR CHANNELS ASSIGNED 


Fig. 8—Channel usage on reuse ring for two-dimensional ring strategy. 


as the percentage of the total number of channel assignments made. 
For example, at a 6.65 percent blocking rate, the channel assigned was 
in use in three coverage areas on the reuse ring for about 30 percent 
of the total channel assignments made. It is interesting to note that the 
assigned channel was being used in the maximum possible number (8) 
of coverage areas on the reuse ring in ine than 0.1 percent of the 
assignments made. 


4.2 Orthogonal NN Strategy 


Some characteristics of the behavior of the two-dimensional ortho- 
gonal NN channel assignment strategy are shown in Fig. 9 as a function 
of the blocking rate of new call-attempts. The curve through the open 
squares shows that for very few of the channel assignments attempted 
there was no usage of the channel assigned in any coverage area on 
the reuse ring. The curve through the open triangles indicates that 
channel assignment was made on the basis of the maximum usage on 
only one side of the reuse ring for 40 to 50 percent of the assignments 
attempted. Below a 10 percent blocking rate, about 50 percent of the 
channel assignments attempted benefitted from orthogonal side optimi- 
zation as indicated by the curves through the open circles and the X’s. 
The curve through the X’s indicates that up to about 25 percent of the 
assignments attempted could have benefitted from usage maximization 
on the third side. 


1628 THE BELL SYSTEM TECHNICAL JOURNAL, SEPTEMBER 1972 


EXPONENTIAL CALL DURATION DISTRIBUTION 


ONLY ONE CHANNEL WITH SIDE 
7 WITH MAXIMUM USAGE 


ONLY ONE CHANNEL WITH ORTHOGONAL 
7 SIDE WITH MAXIMUM USAGE 


MULTIPLE CHANNELS WITH MAXIMUM 
77 USAGE ON ORTHOGONAL SIDE 


/NO USAGE ON REUSE RING 
NO USAGE ON 


7 ORTHOGONAL SIDE 
Se. 


PER CENT OF TOTAL ASSIGNMENTS ATTEMPTED 





1 2 4 6 8 10 20 40 60 80 100 
CALLS BLOCKED PERCENT OF NEW CALL ATTEMPTS 


Fig. 9—Performance characteristics of two-dimensional NN channel assignment. 


Vv. CONCLUSIONS 


The results from the simulation of the two-dimensional dynamic 
channel assignment mobile radio systems show that these systems 
operate at very low blocking until the traffic offered reaches some 
critical value. Small increases in loading above this value produce a 
considerable increase in the blocking of new call attempts and result 
in very little increase in the traffic carried by the system. The loading 
at which blocking begins to occur in the two-dimensional systems is 
somewhat greater than the loading at which blocking begins to occur 
in one-dimensional dynamic channel assignment systems, and is con- 
siderably greater than the loading at which a fixed channel assignment 
system begins to incur significant blocking. For example, at a blocking 
rate, B of 1 percent, the traffic carried, TC, in Erlangs per channel per 
coverage area and the traffic offered, TO, in Erlangs per channel per 
coverage area (TO = TC/(1i — B)) for the systems studied which had an 
average of 10 channels available per coverage area and a reuse interval 
of 4 are: two-dimensional dynamic channel assignment TO = 0.636, 
TC = 0.630; one-dimensional dynamic channel assignment systems, 
TO = 0.631, TC = 0.625; fixed channel assignment systems, TO = 
0.444, TC = 0.440. 

It was found that the two strategies aimed at channel reuse optimiza- 
tion performed similarly and that both were somewhat better (carrying 
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approximately 5 percent more traffic at a given blocking rate) than a 
channel assignment strategy which chose the first available channel 
encountered in a channel search without any regard to reuse optimiza- 
tion. The reason that these dynamic channel assignment systems do not 
carry as much traffic at high blocking rates as fixed channel assignment 
systems is that they are not able to maximize channel reuse as they 
serve the randomly-offered call attempts. It should be possible to re- 
assign channels to calls in progress so that the channel reuse is more 
concentrated and thus more traffic could be carried per channel at a 
given blocking rate. 
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The Buried Channel Charge Coupled Device 


By R. H. WALDEN, R. H. KRAMBECK, R. J. STRAIN, 
J. McKENNA, N. L. SCHRYER, and G. E. SMITH 


(Manuscript received May 18, 1972) 


The charge-coupled device, as described originally by Boyle and 
Smith,’ operates by moving minority carriers along the surface of a 
semiconductor with voltage pulses applied to metal electrodes which 
are separated from the semiconductor by an insulating layer. The transit 
from one electrode to the next is determined by the minority carrier 
transport under the influence of their own potential, fringing fields, 
and diffusion, and by the trapping properties of interface states. The 
transport limitations are largely determined by device geometry;” for 
10-um electrodes, thermal diffusion is predominantly responsible for 
transferring the last small amounts of charge forward and limits efficient 
operation to clock frequencies below 10 MHz. Surface state trapping is 
much less dispersive, and even at low frequencies 10"’ states/em’” eV can 
impose the requirement for regeneration after as few as 100 transfers.*"* 
In order to circumvent these problems, Boyle and Smith’ have proposed 
a modified CCD structure in which the charges do not flow at the semi- 
conductor surface; instead they are confined to a channel which lies 
beneath the surface. This buried channel device has the potential of 
eliminating surface state trapping. (Bulk trapping should be several 
orders of magnitude less important as a CCD loss mechanism.) Calcula- 
tions show that this modification will give rise to increased fringing 
fields under the CCD electrodes. Diffusion is replaced as an important 
factor in the intrinsic transfer process by the more powerful field-aided 
transfer. This leads to fast, efficient transport even when very little 
charge remains to be transferred. In addition to these advantages, the 
higher mobility found in the bulk of the semiconductor should further 
enhance the speed of the device.” 

The proposed buried channel structure is shown in Fig. 1 as a three- 
phase, p-channel CCD. It consists of a p-type layer of Si on an n-type 
Si substrate with p* contacts at either end of the channel. The structure 
is completed with an SiO, film and metal electrodes. The thickness of 
the p-layer is a few um and the acceptor concentration is 10-20 times 
the donor density of the substrate. In order to form a buried channel, 
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Fig. 1—Cross-sectional view of buried channel CCD. 


it is necessary to completely deplete the p-type layer and part of the 
adjacent n-region. This is accomplished by putting the output diode 
at a voltage V, , a reverse bias strong enough to drain all of the holes out 
of the p-layer (a total charge of the order of 10'* cm™*). A solution of 
the one-dimensional Poisson’s equation for this depleted condition is 
shown in Fig. 2, where it was assumed that the p-layer has a uniform 
concentration of 2 X 10'° cm’’ and is 5 um thick, and the SiO, is 
1000 A thick. The substrate is taken to be uniformly doped. to 
1 X 10'* em™*. The graph shows negative potential plotted upward as 
a function of distance away from the Si-SiO, interface for an applied 
gate voltage of 0 V. The peak at —37 V in the curve represents a poten- 
tial minimum for holes, and its position, about 4 um from the interface, 
is the buried channel. This value of peak voltage also gives some idea 
of the potential that must be applied to the output diode in order to 
deplete the p-layer. 

The charge-carrying capability of this structure has been estimated 
to be about 4 X 107° C/cm’ and the effect of this charge on the potential 
is shown in Tig. 2. A surface device would accommodate a charge load 
of about 3 X 107’ C/cm’, almost an order of magnitude more, and this 
is one of the penalties of using a buried channel structure, a reduced 
charge capacity. 

In order to ascertain the fringing fields in the channel, two-dimensional 
potential calculations’ have been made numnierically using the buried 
channel structure. The first structure considered for this calculation, 
corresponding to the experimental device described below, is shown in 
Fig. 3 and consisted of 45-um-wide plates with 5-um interelectrode gaps. 
The channel potential curve in Fig: 3 shows a rather disturbing feature, 
i.e., the existence of potential wells under the interelectrode spaces. 
It is understandable that such wells will exist since the floating region 
between the plates would go to —800 V without the influence of fringing 
fields from the grounded plates. These wells will, of course, fill with 
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Fig. 2—Potential as a function of distance into the silicon showing an empty and a 
filled well. 


mobile charge, but the amount contained in each is going to depend on 
the potentials of the spaces under the neighboring CCD plates. This 
means that the wells can exchange charge with any signal that might 
be present, leading to a severe loss. It would be desirable then to have 
the interelectrode potential vary monotonically across the gap; one way 
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Fig. 3—Two-dimensional calculation of potential along the buried channel. The 
finite interelectrode gaps give rise to the potential minima. 
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Fig. 4—Two-dimensional calculation of potential along the buried channel for 
the case of zero spacing. Also shown is potential plot for a surface device for com- 
parison. 


to do this would be to let the spacing go to zero. The potential calcula- 
tion was repeated for the configuration shown in Fig. 4 for a plate 
length of 10 wm and for a zero spacing between metals. Two curves are 
shown, one for the buried channel CCD and one for a surface device 
for comparison. Here, the undesirable charge pockets are clearly missing 
and there is a significant amount of tilt in the potential under the plate 
at —10 V (the one where charge would be transferred from). The surface 
device, on the other hand, shows a relatively flat potential. The field 
under the —10-V plate for the buried channel is at least 3 X 10° V/cm, 
while the minimum field for the surface device is about 10 times smaller. 
Transport calculations, using this field value and the indicated plate 
length, show that only 107’ of the charge will remain behind after 1 ns. 
Thus operation of the shift register at a clock frequency of several 
hundred MHz should be possible. 

In order to examine these ideas experimentally, some test structures 
were fabricated. The p-layer was formed by a shallow ion implantation 
_ of 10” boron atoms/cm’ followed by a 7-hour diffusion at 1200°C. It 
was estimated that the final junction depth was Sia 7 um. The substrate 
had a doping density of 1 X 10'* cm™* and 1000 A of SiO, was employed 
as the insulator. Aluminum metallizations were defined using standard 
photoresist techniques to yield the CCD plate structure shown in Fig. 3, 
i.e., 45-um plates and 5-um spaces. These devices failed to operate as 
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CCD’s because of the gap problem which was subsequently understood 
as a result of the calculation. However, the existence of a buried channel 
was verified by operating the device as an IGFET. This test was made 
by tying all of the electrodes to a single gate voltage. The channel 
conductance was then measured as a function of the source and drain 
voltages. This experiment was repeated at a series of different gate 
voltages and the resulting cutoff voltages are plotted against gate 
voltage in Fig. 5. The solid curve was obtained by means of the same 
potential calculation that was used for Fig. 3 and the agreement between 
theory and experiment is good. 

In conelusion, a buried channel CCD has been described which offers 
the advantages of improved high-frequency response relative to surface 
devices because of stronger fringe-field-aided transfer and also because 
carrier mobilities have higher values in bulk Si than at the Si-Si0, 
interface. Further, the transfer efficiency should also be better at lower 
frequencies since surface state losses are eliminated. One significant 
design criterion has been established: the interelectrode potential in the 
channel of a polyphase device must vary monotonically across the gap. 
This can be achieved by making the gap very narrow, e.g., by using an 
overlapping metallization structure’’*® or by undercut isolation tech- 
niques.” Once the gap problem is alleviated, the major remaining 
disadvantage of the buried channel device is its reduced charge-carrying 
capability. It is possible, however, to trade some degree of fringe-field- 
aided transfer for charge capacity by moving the channel closer to the 
semiconductor surface. 
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Fig. 5—Experimental plot of channel cutoff voltage as a function of gate voltage. 
Also shown is a theoretical curve for an average doping density in the p-layer of 
1.4 X 104 cm™. 
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